R code for Multiple Isotonic Regression
# MULTIPLE ISOTONIC REGRESSION

To fit a data set (**x**_{i},y_{i}), i=1,...,n, where **x**_{i}
are in R^{p}. It is assumed that the expected value of y is increasing in **x**.
This means that if m=E(y), we have m_{i}<= m_{j} if
**x**_{i} <=**x**_{j} where the last inequality is dimension-wise.
Inputs to the code are:

xmatrix: a matrix of predictors with n rows and p columns;

ys: a vector of responses;

xpred: a matrix with p columns (distinct points);

decr: which should be FALSE if increasing isotonic regression is desired and TRUE if
decreasing isotonic regression is desired.
Outputs are:

theta: estimated response;

thetap: estimated response at rows of xpred

df: degrees of freedom of the fit

### R code for multiple isotonic regression

If p=2, then the function persp may be used to plot the fit. Here is some example code:
x1=runif(100);x2=runif(100)

xmatrix=cbind(x1,x2)

ys=x1+exp(x2)+rnorm(100)/4

x1max=max(x1);x2max=max(x2)

x1min=min(x1);x2min=min(x2)

xp1=(trunc(0:399/20))/19*(x1max-x1min)+x1min

xp2=(0:399-19*20*xp1)/19*(x2max-x2min)+x2min

xpred=cbind(xp1,xp2)

x1g=0:19/19;x2g=x1g

ans=multiso(xmatrix,ys,xpred,FALSE)

surf=matrix(ans$thetap,nrow=20,ncol=20)

persp(x1g,x2g,surf,theta=-30,tick="detailed")

To test constant versus isotonic regression function, use the following code:
### R code for test of constant vs isotonic regression

inputs are as above, with the addition of nloop, which is the number of simulations
used to get the mixing distribution for the test. If n is in the hundreds, this might take
a while for nloop =1000 or larger.
Outputs are as above with the addition of pval, to test null =constant.

To test constant versus isotonic regression function in the presence of a categorical covariate,
use the following code:
### R code for test of constant vs isotonic regression with categorical covariate

inputs are as above, with the addition of zvec, the vector of values of the covariate, which must
range from 1 to z.
Outputs are as above with the addition of pval, to test null =constant.

meyer@stat.colostate.edu