Code Catalog

A one-page catalog of working scripts and quickstarts for np, npRmpi, and crs.
Keywords

np, npRmpi, crs, code examples, scripts, quickstart

This page is meant to be the quickest route to working scripts. It groups the existing code by package, task, and runtime mode so that you can scan, copy, and then drill down only if you need more context.

How to use this page

  • If you want only the smallest starter scripts, go first to Quickstarts.
  • If you want something you can paste into R immediately, start with the short blocks below.
  • If you want a fuller script, use the tables further down.
  • If you know a function name, use the site search box or the Reference page.

Source and scope conventions

The catalog combines Gallery-authored quickstarts and simulations, examples based on documented package datasets, and adaptations of published research. The substantive research examples identify the original paper and reference software, state which calculations are supplied by np, and disclose material departures from the source implementation. Short teaching examples cite their data source on the corresponding page; they are not presented as reproductions of published estimates.

Minimal starters you can copy immediately

np: kernel regression

Source file: np_regression_quickstart.R

## Minimal np regression example.
##
## The intended workflow is:
## 1. compute a bandwidth object,
## 2. fit the regression estimator,
## 3. inspect the result and a simple fitted curve.

library(np)
options(np.messages = FALSE)

data(cps71, package = "np")
dat <- cps71[, c("logwage", "age")]

## Select the bandwidth object, then fit the regression.
bw <- npregbw(logwage ~ age, data = dat, regtype = "ll", bwmethod = "cv.aic")
fit <- npreg(bws = bw, data = dat)

## Review the bandwidth choice and the fitted regression.
summary(bw)
summary(fit)

## Let the plot method handle the basic fitted-curve display.
plot(fit)

Start here if you want the core np workflow in its simplest form. For more context, see Kernel Primer or Quickstarts.

np: density estimation

Source file: np_density_quickstart.R

## Minimal np density-estimation example.
##
## The intended workflow is:
## 1. compute a bandwidth object,
## 2. fit the density estimator,
## 3. inspect the result.

library(np)
options(np.messages = FALSE)

## Start with one variable so the density workflow is easy to read.
data(faithful, package = "datasets")
dat <- data.frame(waiting = faithful$waiting)

## Select the smoothing parameter, then fit the density.
bw <- npudensbw(~ waiting, data = dat, bwmethod = "cv.ml")
fhat <- npudens(bws = bw, data = dat)

## Review the bandwidth choice and the fitted density object.
summary(bw)
summary(fhat)

Start here if the object of interest is a density rather than a regression mean. For distribution, conditional-density, conditional-distribution, and quantile starters, see Quickstarts and Density, Distribution, Quantiles.

npRmpi: same workflow, MPI initialized once

Source file: nprmpi_session_quickstart.R

## Minimal modern npRmpi example.
##
## The intended workflow is:
## 1. initialize MPI once in session/spawn mode,
## 2. write ordinary np-style code,
## 3. quit cleanly at the end.

library(npRmpi)

## Initialize once, then work as if this were an ordinary np script.
npRmpi.init(nslaves = 1)

set.seed(1)
x <- runif(200)
y <- sin(2 * pi * x) + rnorm(200, sd = 0.2)
dat <- data.frame(y, x)

## Fit the bandwidth and regression objects.
bw <- npregbw(y ~ x, regtype = "ll", bwmethod = "cv.ls", data = dat)
fit <- npreg(bws = bw, data = dat)

## Review the objects, plot the fit, then exit cleanly.
summary(bw)
summary(fit)

plot(fit)

npRmpi.quit()

This is the modern default instructional path for npRmpi: initialize MPI once, run ordinary np-style code, then quit cleanly. For mode details, see MPI and Large Data.

crs: spline regression

Source file: crs_quickstart.R

## Minimal crs spline-regression example.
##
## This is the smallest useful workflow:
## 1. fit a spline model,
## 2. inspect the summary,
## 3. optionally move on to plotting or tighter search control.

library(crs)
options(crs.messages = FALSE)

## Generate a simple nonlinear surface with two predictors.
set.seed(42)
n <- 250
x1 <- runif(n)
x2 <- runif(n)
y <- sin(2 * pi * x1) + x2 + rnorm(n, sd = 0.2)
dat <- data.frame(y, x1, x2)

## Fit the spline model, inspect it, then let plot.crs() build the surface.
fit <- crs(y ~ x1 + x2, data = dat)

summary(fit)
plot(fit, perspective = TRUE, view = "fixed")

Use this when regression splines are the natural first tool. For the conceptual route, see Spline Primer.

np: kernel methods in serial

Script Main topic Functions/features Start reading here
np_regression_quickstart.R Minimal kernel-regression workflow npregbw, npreg Kernel Primer
np_regression_nomad_quickstart.R Minimal LP automatic-search workflow npreg, nomad = "auto" Kernel Primer
np_density_quickstart.R Minimal density-estimation workflow npudensbw, npudens Density, Distribution, Quantiles
np_density_bounds_quickstart.R Known and unknown density support bounds npudensbw, npudens, npcdensbw, npcdens Choosing support bounds
np_distribution_quickstart.R Minimal distribution-estimation workflow npudistbw, npudist Density, Distribution, Quantiles
np_conditional_density_quickstart.R Minimal conditional-density workflow npcdensbw, npcdens Density, Distribution, Quantiles
np_conditional_distribution_quickstart.R Minimal conditional-distribution workflow npcdistbw, npcdist Density, Distribution, Quantiles
np_quantile_quickstart.R Minimal quantile-regression workflow npcdistbw, npqreg Density, Distribution, Quantiles
np_lsqreg_quickstart.R Minimal location-scale quantile workflow nplsqregbw, nplsqreg Density, Distribution, Quantiles
np_copula_quickstart.R Minimal copula-estimation workflow npudistbw, npcopula Density, Distribution, Quantiles
np_classification_quickstart.R Minimal classification / conditional-mode workflow npconmode Classification and Modes
np_plotting_quickstart.R Minimal plotting and interval workflow npreg, plot, predict Plotting and Intervals
np_entropy_quickstart.R Minimal entropy/testing workflow npunitest Entropy and Testing
np_semiparametric_quickstart.R Minimal semiparametric workflow npplreg Semiparametric Models
np_instrumental_regression_quickstart.R Minimal instrumental-regression workflow npregiv Kernel Methods
np_significance_quickstart.R Minimal significance-testing workflow npreg, npsigtest Significance and Specification
np_specification_quickstart.R Minimal model-specification workflow lm, npcmstest Significance and Specification

Additional serial scripts

Script Main topic Functions/features Start reading here
regression_intro_a.R Univariate regression and plotting npreg, local constant, local linear Kernel Primer
regression_intro_b.R Derivative estimation npreg, derivative output Plotting and Intervals
regression_intro_c.R Plotting fitted objects and intervals plot(), asymptotic and bootstrap error handling Plotting and Intervals
regression_multivar_a.R Multivariate regression mixed predictors, multivariate fit Multivariate Regression and Prediction
sharp_rd_inference.R Sharp regression discontinuity npreg, RD bandwidths, bias-corrected inference Causal Inference
cholesterol_simultaneous.R Randomized treatment effects npreg, simultaneous grid bands, joint test Causal Inference
np_crossfit_cate.R Doubly robust treatment effects npreg, lasso nuisances, cross-fitting, multiplier bands Doubly Robust Treatment Effects
np_composition_did.R Difference-in-differences with changing composition npconmode, npreg, fresh nuisance estimation, clustered ATT inference Difference-in-Differences

np: practitioner-oriented extras

These are compact workflows drawn from the practitioner-oriented material in the Cambridge manuscript. They are useful when the applied task is not just regression fitting.

Script Main topic Functions/features Start reading here
np_probability_quickstart.R Smooth probability estimation on discrete support npudens, Aitchison-Aitken kernel Density, Distribution, Quantiles
np_density_resample_quickstart.R Smooth resampling from a fitted density npudensbw, npudens, Gaussian-kernel resampling Density, Distribution, Quantiles
np_conditional_density_resample_quickstart.R Smooth resampling from a fitted conditional density npcdensbw, npcdens, weighted conditional resampling Density, Distribution, Quantiles

npRmpi: modern session-mode quick start

If you want one current downloadable script before anything else, start here.

Script Main topic Functions/features Start reading here
nprmpi_session_quickstart.R Modern session / spawn workflow npRmpi.init, npregbw, npreg MPI and Large Data
nprmpi_session_nomad_quickstart.R Modern session LP automatic-search workflow npRmpi.init, npreg, nomad = "auto" MPI and Large Data
nprmpi_attach_quickstart.R Modern attach-mode workflow npRmpi.init(mode = \attach\), npregbw, npreg MPI and Large Data
nprmpi_profile_quickstart.R Modern profile/manual-broadcast workflow mpi.bcast.cmd, np.mpi.initialize, npregbw, npreg MPI and Large Data

Inline source for common scripts

The blocks below are pulled directly from the script files in www/, so they stay aligned with the downloadable versions.

Source file: np_regression_quickstart.R

## Minimal np regression example.
##
## The intended workflow is:
## 1. compute a bandwidth object,
## 2. fit the regression estimator,
## 3. inspect the result and a simple fitted curve.

library(np)
options(np.messages = FALSE)

data(cps71, package = "np")
dat <- cps71[, c("logwage", "age")]

## Select the bandwidth object, then fit the regression.
bw <- npregbw(logwage ~ age, data = dat, regtype = "ll", bwmethod = "cv.aic")
fit <- npreg(bws = bw, data = dat)

## Review the bandwidth choice and the fitted regression.
summary(bw)
summary(fit)

## Let the plot method handle the basic fitted-curve display.
plot(fit)

Source file: np_regression_nomad_quickstart.R

library(np)
data(cps71, package = "np")

## Use the shortest modern LP route.
## nomad = "auto" uses exhaustive search for p = 1 and NOMAD otherwise.
fit_lp <- npreg(logwage ~ age, data = cps71, nomad = "auto")

## Plot first, then inspect the fitted object in text form.
plot(fit_lp)
summary(fit_lp)

Source file: np_density_quickstart.R

## Minimal np density-estimation example.
##
## The intended workflow is:
## 1. compute a bandwidth object,
## 2. fit the density estimator,
## 3. inspect the result.

library(np)
options(np.messages = FALSE)

## Start with one variable so the density workflow is easy to read.
data(faithful, package = "datasets")
dat <- data.frame(waiting = faithful$waiting)

## Select the smoothing parameter, then fit the density.
bw <- npudensbw(~ waiting, data = dat, bwmethod = "cv.ml")
fhat <- npudens(bws = bw, data = dat)

## Review the bandwidth choice and the fitted density object.
summary(bw)
summary(fhat)

Source file: np_density_bounds_quickstart.R

## Known support bounds with the default Gaussian kernel (np 0.80-1).
## Source this file to fit and plot all six seed-42 examples.
## Save the plot with:
## png("density_bounds.png", width = 1800, height = 1400, res = 180)
## plot_density_bounds(bounds_examples, main = "One sample | n = 250 | seed = 42")
## dev.off()
## The separate 100-draw Gallery averages are not run by this quickstart.

library(np)
options(np.messages = FALSE)
RNGkind("Mersenne-Twister", "Inversion", "Rejection")
n <- 250L
seed <- 42L
search_seed <- 100000L + seed
designs <- c("proportion", "duration", "gaussian")
bounds_examples <- list()

for (design in designs) {
  bounds <- switch(design, proportion = c(0, 1),
                   duration = c(0, Inf), gaussian = c(-Inf, Inf))
  yy <- switch(design, proportion = seq(-.2, 1.2, length.out = 281),
               duration = seq(-.8, 8, length.out = 221),
               gaussian = seq(-4.5, 4.5, length.out = 201))
  inside <- yy >= bounds[1] & yy <= bounds[2]

  ## Unconditional: the same sample, separate CVML bandwidth searches.
  set.seed(seed)
  y <- switch(design, proportion = runif(n),
              duration = rexp(n), gaussian = rnorm(n))
  set.seed(search_seed); npseed(search_seed)
  b0 <- npudensbw(dat = data.frame(y), bwmethod = "cv.ml", nmulti = 5)
  set.seed(search_seed); npseed(search_seed)
  b1 <- npudensbw(dat = data.frame(y), bwmethod = "cv.ml", nmulti = 5,
                  ckerbound = "fixed", ckerlb = bounds[1], ckerub = bounds[2])
  u0 <- fitted(npudens(bws = b0, tdat = data.frame(y), edat = data.frame(y = yy)))
  ## Public evaluation is inside the declared support. Zero extension is
  ## an explicit plotting convention; no in-support estimates are changed.
  u1 <- numeric(length(yy))
  u1[inside] <- fitted(npudens(bws = b1, tdat = data.frame(y),
                             edat = data.frame(y = yy[inside])))
  tu <- switch(design, proportion = dunif(yy),
               duration = dexp(yy), gaussian = dnorm(yy))

  ## Conditional: one continuous X and one continuous Y, default LC.
  set.seed(seed)
  x <- switch(design, proportion = runif(n),
              duration = rexp(n), gaussian = rnorm(n))
  y <- switch(design,
    proportion = {a <- .8 * (2*x - 1); u <- runif(n)
                  2*u / ((1-a) + sqrt((1-a)^2 + 4*a*u))},
    duration = rexp(n, rate = 1 / (1 + .5*x)),
    gaussian = rnorm(n, mean = x))
  ## Display the conditional density at the population median of X.
  x0 <- switch(design, proportion = .5, duration = qexp(.5), gaussian = 0)
  grid <- data.frame(x = rep(x0, length(yy)), y = yy)
  set.seed(search_seed); npseed(search_seed)
  b0c <- npcdensbw(xdat = data.frame(x), ydat = data.frame(y),
                   bwmethod = "cv.ml", nmulti = 5)
  set.seed(search_seed); npseed(search_seed)
  b1c <- npcdensbw(xdat = data.frame(x), ydat = data.frame(y),
                   bwmethod = "cv.ml", nmulti = 5,
                   cxkerbound = "fixed", cxkerlb = bounds[1], cxkerub = bounds[2],
                   cykerbound = "fixed", cykerlb = bounds[1], cykerub = bounds[2])
  c0 <- fitted(npcdens(bws = b0c, txdat = data.frame(x), tydat = data.frame(y),
                      exdat = grid["x"], eydat = grid["y"], proper = FALSE))
  c1 <- numeric(length(yy))
  c1[inside] <- fitted(npcdens(bws = b1c, txdat = data.frame(x), tydat = data.frame(y),
                              exdat = grid[inside, "x", drop = FALSE],
                              eydat = grid[inside, "y", drop = FALSE], proper = FALSE))
  tc <- switch(design,
    proportion = (1 + .8*(2*x0 - 1)*(2*yy - 1)) * inside,
    duration = dexp(yy, rate = 1 / (1 + .5*x0)), gaussian = dnorm(yy, mean = x0))

  if (design == "gaussian") {
    ## Test the estimator separately from bandwidth-search rounding.
    same <- npcdens(bws = c(b0c$ybw, b0c$xbw),
                    txdat = data.frame(x), tydat = data.frame(y),
                    exdat = grid["x"], eydat = grid["y"], proper = FALSE,
                    cxkerbound = "fixed", cxkerlb = -Inf, cxkerub = Inf,
                    cykerbound = "fixed", cykerlb = -Inf, cykerub = Inf)
    stopifnot(identical(u0, u1), identical(c0, fitted(same)))
    print(c(unconditional_difference = max(abs(u0 - u1)),
            conditional_same_bandwidth = max(abs(c0 - fitted(same))),
            conditional_separate_cvml = max(abs(c0 - c1))))
  }
  bounds_examples[[design]] <- list(y = yy, bounds = bounds, x0 = x0,
    unconditional = data.frame(truth = tu, ordinary = u0, fixed = u1),
    conditional = data.frame(truth = tc, ordinary = c0, fixed = c1))
}

plot_density_bounds <- function(examples, main = "One sample",
                                subtitle = "One fitted estimate per curve",
                                ylimits = NULL) {
  op <- par(no.readonly = TRUE)
  on.exit(par(op))
  par(mfrow = c(2, 3), cex = 1, mar = c(4, 4, 3.5, .8), oma = c(0, 0, 6.5, 0))
  cols <- c("#222222", "#D55E00", "#0072B2")
  titles <- c(proportion = "Proportion: two endpoints",
              duration = "Duration: one endpoint", gaussian = "Gaussian: no finite endpoint")
  for (family in c("unconditional", "conditional")) for (design in names(examples)) {
    z <- examples[[design]]; v <- z[[family]]
    inside <- z$y >= z$bounds[1] & z$y <= z$bounds[2]
    ylim <- if (is.null(ylimits)) c(0, max(v)*1.07) else ylimits[[family]][[design]]
    plot(z$y, v$truth, type = "n", ylim = ylim,
         xlab = "Value of Y", ylab = "Density", main = titles[design], cex.main = .9)
    mtext(if (family == "unconditional") "Unconditional" else
            sprintf("Conditional at X = %.2f", z$x0), side = 3, line = .3, cex = .85)
    if (is.finite(z$bounds[1])) rect(min(z$y), 0, z$bounds[1], par("usr")[4],
                                   col = "#F0F0F0", border = NA)
    if (is.finite(z$bounds[2])) rect(z$bounds[2], 0, max(z$y), par("usr")[4],
                                   col = "#F0F0F0", border = NA)
    ## Break at support endpoints rather than draw a spurious sloping join.
    for (ix in split(seq_along(z$y), cumsum(c(TRUE, diff(inside) != 0)))) {
      lines(z$y[ix], v$truth[ix], col = cols[1], lwd = 2.4)
      lines(z$y[ix], v$fixed[ix], col = cols[3], lwd = 2.4, lty = 3)
    }
    lines(z$y, v$ordinary, col = cols[2], lwd = 2, lty = 2)
    if (design == "gaussian") lines(z$y, v$fixed, col = cols[3], lwd = 2.4, lty = 3)
  }
  mtext(main, side = 3, outer = TRUE, line = 3.5, font = 2, cex = 1.25)
  mtext(subtitle, side = 3, outer = TRUE, line = 2, cex = .9)
  top <- par("omi")[3] / par("din")[2]
  par(fig = c(0, 1, 1-top, 1), new = TRUE, mar = rep(0, 4), oma = rep(0, 4))
  plot.new()
  legend("bottom", c("True density", "Ordinary Gaussian kernel", "Gaussian kernel with fixed support"),
         col = cols, lty = c(1, 2, 3), lwd = 2.4, bty = "n", horiz = TRUE, cex = .9)
}

plot_density_bounds(bounds_examples, main = sprintf("One sample | n = %d | seed = %d", n, seed))

Source file: np_distribution_quickstart.R

## Minimal np distribution example.
##
## This mirrors the usual density workflow but targets the
## unconditional distribution function instead.

library(np)
options(np.messages = FALSE)

## Mirror the basic density setup, but target the distribution instead.
data(faithful, package = "datasets")
dat <- data.frame(waiting = faithful$waiting)

## Select the bandwidth, then fit the unconditional distribution.
bw <- npudistbw(~ waiting, data = dat, nmulti = 1)
Fhat <- npudist(bws = bw, data = dat)

## Review the smoothing choice and the fitted distribution object.
summary(bw)
summary(Fhat)

Source file: np_conditional_density_quickstart.R

## Minimal np conditional-density example.
##
## This keeps the first run small enough to be practical while still
## showing the standard two-step workflow:
## 1. compute a bandwidth object,
## 2. fit the conditional-density estimator.

library(np)
options(np.messages = FALSE)

## Keep the first run small enough that the two-step workflow stays quick.
data(faithful, package = "datasets")
dat <- faithful[seq_len(120), c("eruptions", "waiting")]

## Select the bandwidth first, then fit the conditional density.
bw <- npcdensbw(eruptions ~ waiting, data = dat, nmulti = 1)
fhat <- npcdens(bws = bw, data = dat)

## Review both the smoothing choice and the fitted object.
summary(bw)
summary(fhat)

Source file: np_conditional_distribution_quickstart.R

## Minimal np conditional-distribution example.
##
## This uses a small first run so the standard two-step workflow stays
## copyable and practical.

library(np)
options(np.messages = FALSE)

## Work on a smaller slice so the first run stays easy to inspect.
data(faithful, package = "datasets")
dat <- faithful[seq_len(120), c("eruptions", "waiting")]

## Select the bandwidth, then fit the conditional distribution.
bw <- npcdistbw(eruptions ~ waiting, data = dat, nmulti = 1)
Fhat <- npcdist(bws = bw, data = dat)

## Review the bandwidth object and the fitted distribution.
summary(bw)
summary(Fhat)

Source file: np_quantile_quickstart.R

## Minimal np quantile-regression example.
##
## The key idea is to compute/select the conditional distribution once
## and then extract one or more quantiles from it.

library(np)
options(np.messages = FALSE)

## Keep the first quantile run compact and reproducible.
data(faithful, package = "datasets")
dat <- faithful[seq_len(120), c("eruptions", "waiting")]

## Build one conditional-distribution fit and request several quantiles.
qfit <- npqreg(eruptions ~ waiting,
               data = dat,
               tau = c(0.25, 0.50, 0.75),
               nmulti = 1)

## Review the bandwidth choice and one representative fitted quantile.
summary(qfit$bws)
summary(qfit)

## Plot the overlaid fitted quantile curves.
plot(qfit)

Source file: np_lsqreg_quickstart.R

## Minimal np location-scale quantile example.
##
## Fit several location-scale quantiles in one call, then inspect the
## selected bandwidth state and compare the fitted curves.

library(np)
options(np.messages = FALSE)

## Keep the first location-scale run compact and reproducible.
data(Italy, package = "np")

## Fit several conditional quantiles using the location-scale route.
qfit <- nplsqreg(gdp ~ ordered(year),
                 data = Italy,
                 tau = c(0.25, 0.50, 0.75),
                 nmulti = 1)

## Review the selected search state and fitted object.
summary(qfit$bws)
summary(qfit)

## Predict the stored median without repeating bandwidth selection.
qmed <- predict(qfit, tau = 0.50)
head(qmed)

## Plot the fitted location-scale quantile curves.
plot(qfit)

Source file: np_copula_quickstart.R

## Minimal np copula example.
##
## This uses the direct formula interface and lets npcopula() build
## the corresponding marginal distribution bandwidth object internally.

library(np)
library(MASS)

set.seed(42)
n <- 1000
rho <- 0.95
Sigma <- matrix(c(1, rho, rho, 1), 2, 2)

## Simulate a simple bivariate sample with strong dependence.
dat <- as.data.frame(mvrnorm(n = n, mu = c(0, 0), Sigma = Sigma))
names(dat) <- c("x", "y")

## Let npcopula() create a small default probability grid.
copula_fit <- npcopula(~ x + y, data = dat, neval = 20)

## Inspect the fitted copula object and the retained bandwidth object.
summary(copula_fit)
summary(copula_fit$bws)

## Evaluate on a small Cartesian grid of marginal probability coordinates.
u_eval <- list(u1 = c(0.25, 0.50, 0.75),
  u2 = c(0.25, 0.50, 0.75))
predict(copula_fit, u = u_eval)

Source file: np_classification_quickstart.R

## Minimal np classification / conditional-mode example.
##
## This keeps the first run compact while still showing the basic
## nonparametric classification route.

library(np)
data(birthwt, package = "MASS")

## Class the discrete fields the way the estimator expects them.
birthwt$low <- factor(birthwt$low)
birthwt$smoke <- factor(birthwt$smoke)
birthwt$race <- factor(birthwt$race)

## Fit the conditional-mode classifier on a small practical formula.
fit <- npconmode(low ~ smoke + race + age + lwt, data = birthwt, nmulti = 1)

## Inspect the fitted object, then the implied confusion matrix.
summary(fit)
fit$confusion.matrix

Source file: np_plotting_quickstart.R

## Minimal np plotting and interval example.
##
## This fits a simple local-linear model and then shows the
## modern bootstrap plotting route with all interval/band types.

library(np)
options(np.messages = FALSE)

data(cps71, package = "np")

## Fit a simple local-linear regression object.
model.ll <- npreg(logwage ~ age, regtype = "ll", data = cps71)

## Save one example plot so the script works in non-interactive sessions too.
plot_path <- file.path(tempdir(), "np_plotting_quickstart.png")
png(plot_path, width = 700, height = 500)
plot(model.ll,
  errors = "bootstrap",
  bootstrap = "inid",
  band = "all")
dev.off()

## Report where the rendered image landed.
cat("Saved plot to:", plot_path, "\n")

Source file: np_entropy_quickstart.R

## Minimal np entropy/testing example.
##
## This uses the univariate density-equality test because it is a
## compact, fast first run.

library(np)
options(np.messages = FALSE)

## Generate two simple samples so the test output is easy to interpret.
set.seed(1234)
n <- 300
x <- rnorm(n)
y <- rnorm(n)

## Run a compact first test without bootstrap overhead.
test_out <- npunitest(x, y, bootstrap = FALSE)

## Inspect the test summary.
summary(test_out)

Source file: np_semiparametric_quickstart.R

## Minimal np semiparametric example.
##
## This uses a partially linear model because it is the clearest
## lightweight entry point into the semiparametric family.

library(np)
options(np.messages = FALSE)

## Build a small partially linear data-generating setup.
set.seed(42)
n <- 200
x1 <- rnorm(n)
z <- runif(n)
y <- 1 + 2 * x1 + sin(2 * pi * z) + rnorm(n, sd = 0.2)
dat <- data.frame(y, x1, z)

## Fit the partially linear model.
fit <- npplreg(y ~ x1 | z, data = dat)

## Inspect the fitted semiparametric object.
summary(fit)

Source file: np_instrumental_regression_quickstart.R

## Minimal np instrumental-regression example.
##
## This mirrors the basic one-endogenous-regressor setup used in the
## package demos while keeping the first run compact.

library(np)

## Build a small one-endogenous-regressor example with one instrument.
set.seed(42)
n <- 80
v <- rnorm(n, sd = 0.27)
eps <- rnorm(n, sd = 0.05)
u <- -0.5 * v + eps
w <- rnorm(n)
z <- 0.2 * w + v
y <- z^2 + u
iv_dat <- data.frame(y, z, w)

## Fit through the 0.80-1 formula/data interface.
fit_iv <- npregiv(y ~ z | w,
  data = iv_dat,
  method = "Tikhonov",
  nmulti = 1)

## Inspect the selected structural fit and its training residuals.
summary(fit_iv)
head(fitted(fit_iv))
head(residuals(fit_iv))

## Use the dedicated IV plot controls for the fitted structural function.
plot(fit_iv, data.overlay = TRUE, data.rug = TRUE)

Source file: np_significance_quickstart.R

## Minimal np significance-testing example.
##
## The idea is to fit a model with one irrelevant regressor and then
## ask whether the nonparametric significance test detects that.

library(np)
options(np.messages = FALSE)

## Generate one relevant and one irrelevant regressor on purpose.
set.seed(42)
n <- 200
z <- factor(rbinom(n, 1, 0.5))
x1 <- rnorm(n)
x2 <- runif(n, -2, 2)
y <- x1 + x2 + rnorm(n)
dat <- data.frame(z, x1, x2, y)

## Fit the model before asking the significance question.
fit <- npreg(y ~ z + x1 + x2,
  regtype = "ll",
  bwmethod = "cv.aic",
  data = dat)

## Use the current B argument for a compact bootstrap run.
test_out <- npsigtest(fit, B = 99)

summary(fit)
summary(test_out)

Source file: np_specification_quickstart.R

## Minimal np model-specification test example.
##
## The idea is to fit a simple linear model to nonlinear data and then
## ask whether the parametric specification looks too restrictive.

library(np)
options(np.messages = FALSE)

## Simulate nonlinear data, then fit an intentionally simple linear model.
set.seed(42)
n <- 120
x <- runif(n, -2, 2)
y <- x + x^2 + rnorm(n, sd = 0.25)

model_ols <- lm(y ~ x, x = TRUE, y = TRUE)
X <- data.frame(x = x)

## Compare the parametric model to the nonparametric alternative.
test_out <- npcmstest(model = model_ols,
  xdat = X,
  ydat = y,
  B = 99,
  nmulti = 1)

## Inspect both the linear fit and the specification test.
summary(model_ols)
summary(test_out)

Source file: regression_intro_a.R

## Compare local-constant and local-linear kernel regressions.
## Data: 1971 Canadian Census sample supplied by Aman Ullah; see
## Pagan and Ullah (1999), Nonparametric Econometrics.

library(np)
data(cps71, package = "np")

## Fit the default local-constant model and display its fitted curve.
model.lc <- npreg(logwage ~ age, data = cps71)
plot(model.lc)

## Fit a local-linear model, then add its curve in sample order.
model.ll <- npreg(logwage ~ age,
  regtype = "ll",
  data = cps71)
age.order <- order(cps71$age)
lines(cps71$age[age.order], fitted(model.ll)[age.order],
  col = 2, lty = 2)

legend("topleft",
  c("Local constant", "Local linear"),
  lty = c(1, 2),
  col = c(1, 2),
  bty = "n")

Source file: regression_multivar_a.R

## Fit and plot a mixed-data multivariate kernel regression.
## Data: 1976 U.S. Current Population Survey sample supplied by Jeffrey M.
## Wooldridge; see Wooldridge (2000), Introductory Econometrics.

library(np)
data(wage1, package = "np")

## Combine one factor with two numeric predictors in a local-linear fit.
model <- npreg(lwage ~ female + educ + exper,
  regtype = "ll",
  data = wage1)

## Compare asymptotic and bootstrap partial-regression views.
plot(model, errors = "asymptotic")
plot(model, errors = "bootstrap")
plot(model, errors = "bootstrap", band = "all")

## Request derivative partial views from the same fitted object.
plot(model, errors = "bootstrap", band = "all", gradients = TRUE)

Source file: nprmpi_session_quickstart.R

## Minimal modern npRmpi example.
##
## The intended workflow is:
## 1. initialize MPI once in session/spawn mode,
## 2. write ordinary np-style code,
## 3. quit cleanly at the end.

library(npRmpi)

## Initialize once, then work as if this were an ordinary np script.
npRmpi.init(nslaves = 1)

set.seed(1)
x <- runif(200)
y <- sin(2 * pi * x) + rnorm(200, sd = 0.2)
dat <- data.frame(y, x)

## Fit the bandwidth and regression objects.
bw <- npregbw(y ~ x, regtype = "ll", bwmethod = "cv.ls", data = dat)
fit <- npreg(bws = bw, data = dat)

## Review the objects, plot the fit, then exit cleanly.
summary(bw)
summary(fit)

plot(fit)

npRmpi.quit()

Source file: nprmpi_session_nomad_quickstart.R

## Minimal modern npRmpi NOMAD example.
##
## The intended workflow is:
## 1. initialize MPI once,
## 2. use the same nomad = "auto" shortcut available in np,
## 3. inspect the normalized shortcut metadata,
## 4. quit cleanly at the end.

if (!requireNamespace("crs", quietly = TRUE)) {
  stop("This quickstart uses nomad = \"auto\", which requires the 'crs' package ",
    "for NOMAD degree search. Install it with install.packages('crs').")
}

library(npRmpi)

## Initialize once, then write ordinary np-style code.
npRmpi.init(nslaves = 1)

set.seed(7)
n <- 120
x <- runif(n, -1, 1)
y <- x + 0.4 * x^2 + rnorm(n, sd = 0.18)
dat <- data.frame(y, x)

## Fit the LP shortcut and inspect the normalized metadata.
## nomad = "auto" uses exhaustive search for p = 1 and NOMAD otherwise.
fit <- npreg(y ~ x, data = dat, nomad = "auto", degree.max = 2L, nmulti = 1L)

fit$bws$nomad.shortcut
summary(fit)

## Plot the fitted object before quitting cleanly.
plot(fit)

npRmpi.quit()

Source file: nprmpi_attach_quickstart.R

## Minimal attach-mode npRmpi example.
##
## Launch with a pre-created MPI world, for example:
##   mpiexec -env R_PROFILE_USER "" -env R_PROFILE "" -n 2 \
##     Rscript --no-save nprmpi_attach_quickstart.R

library(npRmpi)

## Initialize under an already-created MPI world.
npRmpi.init(mode = "attach")

if (mpi.comm.rank(0L) == 0L) {
  ## Build the data and fit on the coordinator rank.
  set.seed(1)
  x <- runif(200)
  y <- sin(2 * pi * x) + rnorm(200, sd = 0.2)
  dat <- data.frame(y, x)

  bw <- npregbw(y ~ x, regtype = "ll", bwmethod = "cv.ls", data = dat)
  fit <- npreg(bws = bw, data = dat)

  ## Inspect the fitted objects, then shut down cleanly.
  summary(bw)
  summary(fit)

  npRmpi.quit(mode = "attach")
}

Source file: nprmpi_profile_quickstart.R

## Minimal profile/manual-broadcast npRmpi example.
##
## Launch with an explicit profile source, for example:
##   RPROFILE=$(Rscript --no-save -e 'cat(system.file("Rprofile", package="npRmpi"))')
##   mpiexec -env R_PROFILE_USER "$RPROFILE" -env R_PROFILE "" -n 2 \
##     Rscript --no-save nprmpi_profile_quickstart.R

## Initialize the workers through the profile route.
invisible(mpi.bcast.cmd(np.mpi.initialize(), caller.execute = TRUE))

## Build the data once on the caller and broadcast it.
set.seed(1)
x <- runif(200)
y <- sin(2 * pi * x) + rnorm(200, sd = 0.2)
dat <- data.frame(y, x)

invisible(mpi.bcast.Robj2slave(dat))

## Fit the bandwidth and regression objects on the workers.
invisible(mpi.bcast.cmd(bw <- npregbw(y ~ x, regtype = "ll", bwmethod = "cv.ls", data = dat),
  caller.execute = TRUE))

invisible(mpi.bcast.cmd(fit <- npreg(bws = bw, data = dat),
  caller.execute = TRUE))

## Inspect the results, then stop the MPI workers.
summary(bw)
summary(fit)

invisible(mpi.bcast.cmd(mpi.quit(), caller.execute = TRUE))

npRmpi: serial and MPI comparison scripts

These are older comparison scripts, but still useful if you want to see side-by-side serial and MPI routes. For new work, start with the session / spawn example above and then come back to these when you want route parity or legacy comparison.

Task Serial MPI Start reading here
Conditional density estimation npcdensls_serial.R npcdensls_npRmpi.R Density, Distribution, Quantiles
Model specification test npcmstest_serial.R npcmstest_npRmpi.R Significance and Specification
Conditional mode estimation npconmode_serial.R npconmode_npRmpi.R Classification and Modes
Density equality test npdeneqtest_serial.R npdeneqtest_npRmpi.R Entropy and Testing
Single index estimation npindexich_serial.R npindexich_npRmpi.R Semiparametric Models
Partially linear regression npplreg_serial.R npplreg_npRmpi.R Semiparametric Models
Local linear AIC bandwidth selection npregllaic_serial.R npregllaic_npRmpi.R Kernel Primer
Serial dependence test npsdeptest_serial.R npsdeptest_npRmpi.R Entropy and Testing
Significance test npsigtest_serial.R npsigtest_npRmpi.R Significance and Specification
Unconditional density estimation npudensml_serial.R npudensml_npRmpi.R Density, Distribution, Quantiles

npRmpi: local polynomial / constrained examples

The older npglpreg route is deprecated. These retained local polynomial scripts use the current npregbw() + npreghat() route instead. For the automatic LP replacement path, use npreg(..., nomad = "auto") in np or npRmpi.

Script Main topic
lp_k1_prodfunc.R Production-function style local polynomial example
lp_k1.R Local polynomial example
lp_radial_mean.R Constrained radial mean example
lp_radial_deriv.R Constrained radial derivative example

crs: spline examples

Current crs plotting follows the same named-control style as np and npRmpi, and the constrained examples use crshat() rather than manual hat-matrix construction.

Script Main topic Notes
crs_quickstart.R Minimal spline-regression workflow plot.crs() surface
radial_rgl.R Bivariate radial spline surface plot(..., renderer = "rgl"); requires rgl
sine_rgl.R Bivariate sine spline surface plot(..., renderer = "rgl"); requires rgl
radial_constrained_mean.R Constrained spline mean estimation crshat() + quadprog
radial_constrained_first_partial.R Constrained first derivative crshat() + quadprog
radial_constrained_second_partial.R Constrained second derivative crshat() + quadprog

Source file: crs_quickstart.R

## Minimal crs spline-regression example.
##
## This is the smallest useful workflow:
## 1. fit a spline model,
## 2. inspect the summary,
## 3. optionally move on to plotting or tighter search control.

library(crs)
options(crs.messages = FALSE)

## Generate a simple nonlinear surface with two predictors.
set.seed(42)
n <- 250
x1 <- runif(n)
x2 <- runif(n)
y <- sin(2 * pi * x1) + x2 + rnorm(n, sd = 0.2)
dat <- data.frame(y, x1, x2)

## Fit the spline model, inspect it, then let plot.crs() build the surface.
fit <- crs(y ~ x1 + x2, data = dat)

summary(fit)
plot(fit, perspective = TRUE, view = "fixed")

Source file: crs_iv_quickstart.R

## Minimal crs instrumental-regression example.
##
## This is a compact first run for spline-based IV estimation.

library(crs)

## Build a small IV example with one endogenous regressor and one instrument.
set.seed(42)
n <- 100
v <- rnorm(n, sd = 0.27)
eps <- rnorm(n, sd = 0.05)
u <- -0.5 * v + eps
w <- rnorm(n)
z <- 0.2 * w + v
y <- z^2 + u
iv_dat <- data.frame(y, z, w)

## Fit through the current formula/data interface.
fit_iv <- crsiv(y ~ z | w,
  data = iv_dat,
  method = "Landweber-Fridman",
  cv = "exhaustive",
  nmulti = 1,
  cv.threshold = 0)

## Inspect the selected structural fit and its training residuals.
summary(fit_iv)
head(fitted(fit_iv))
head(residuals(fit_iv))

Interactive demos for teaching and exploration

Script Main topic Notes
manipulate_density.R Univariate density with sliders RStudio + manipulate
manipulate_eruptions.R Old Faithful density RStudio + manipulate
manipulate_distribution.R Univariate distribution RStudio + manipulate
manipulate_bivariate_density.R Bivariate density RStudio + manipulate
manipulate_faithful_density.R Old Faithful bivariate density RStudio + manipulate
manipulate_bivariate_distribution.R Bivariate distribution RStudio + manipulate
manipulate_wage1.R Multivariate regression partial plots RStudio + manipulate
manipulate_constrained_local_polynomial.R Constrained local polynomial regression requires quadprog
manipulate_constrained_local_polynomial_production.R Constrained production function requires quadprog
manipulate_constrained_spline.R Constrained spline surface requires quadprog
manipulate_constrained_spline_derivative.R Constrained spline derivative requires quadprog
manipulate_copula.R Copula plots requires mnormt
manipulate_socco.R SOCCO nonlinear-system demo requires data.RData and figures.R

Notes

  • The Gallery retains the older teaching scripts, while the current pages route readers to them more directly and supply the missing source context.
  • The short blocks near the top are there so you can copy a minimal working pattern without downloading anything first.
  • For MPI-specific launch semantics, see MPI and Large Data.
  • For short copyable code blocks rather than full scripts, start with Worked Examples.
Back to top