R/logreg.R

logreg <-
function(x,y,xout=FALSE,outfun=outpro,plotit=FALSE,POLY=FALSE,
xlab="X",ylab="Y",zlab="",SCALE=FALSE,expand=.5,theta=50,phi=25,
duplicate="error",ticktype="simple",...){
#
# Perform  logistic regression.
# The predictors are assumed to be stored in the n by p matrix x.
# The y values should be 1 or 0.
#
# xout=T will remove outliers from among the x values and then fit
# the regression line.
#  Default:
# One predictor, a mad-median rule is used.
# With more than one, projection method is used.
#
# outfun=out will use MVE method
#
#  plotit=TRUE will plot regression line
#  POLY=T,  will plot regression line assuming predictor
#  is in  col 1 of x and other columns are x (in col 1) raised to some power
#   or some other function of x
#
x<-as.matrix(x)
p=ncol(x)
xy=elimna(cbind(x,y))
x=xy[,1:ncol(x)]
y=xy[,ncol(xy)]
x<-as.matrix(x)
if(xout){
flag<-outfun(x,...)$keep
x<-x[flag,]
y<-y[flag]
}
x<-as.matrix(x)
if(p==1 || POLY){
xord=order(x[,1])
x=x[xord,]
y=y[xord]
}
fitit=glm(formula=y~x,family=binomial)
init<-summary(fitit)
if(plotit){
vals=fitted.values(fitit)
if(p==1){
plot(x,y,xlab=xlab,ylab=ylab)
lines(x,vals)
}
if(p==2){
if(!scale)print("With dependence, suggest using scale=T")
fitr=vals
iout<-c(1:length(fitr))
nm1<-length(fitr)-1
for(i in 1:nm1){
ip1<-i+1
for(k in ip1:length(fitr))if(sum(x[i,]==x[k,])==2)iout[k]<-0
}
fitr<-fitr[iout>=1] # Eliminate duplicate points in the x-y plane
#                 This is necessary when doing three dimensional plots
#                 with the R function interp
mkeep<-x[iout>=1,]
fit<-interp(mkeep[,1],mkeep[,2],fitr,duplicate=duplicate)
persp(fit,theta=theta,phi=phi,expand=expand,
scale=scale,xlab=xlab,ylab=ylab,zlab=zlab,ticktype=ticktype)
}
}
init$coef
}
musto101/wilcox_R documentation built on May 23, 2019, 10:52 a.m.