Hello, I've written a small function that's supposed to save me some time, and it's ending up killing it- the intention is to iteratively subset a dataset fram on framevec, fit a model (either lm or nls depending on type) and return the r2 or AIC from the model, respectively. Although as far as I can tell in my code the plots are dependent on the fit of the model to the data and the r2 and AIC reported are also dependent on the re-fitted model, the plots show the same linear or non-linear model used for every subset of the data.
However, the r2 and AIC values come back different for each subset. When I do the subsetting and model fitting outside the function, the model fit is different, with different slopes, for each subset of the data. I'm going to end up doing this without the function if I don't solve this soon. Any help much appreciated. #a is the number of times to loop the tenth of percent, b is first dimension of mfrow, frame is the dataframe, framevec is the vector within the dataframe you're using to subset (should be a percentage), form 1 is the simple y~x for the plot, form 2 is y~x for regression, type is lm 1 or 2 nls ,form 3 is the formula for the nls, sx and sn are the start values for nls plotter<-function(a,b,fram,framvec,form1,form2,form3, type=1, xlm, ylm, sx=.01,sn=25){ g<-ceiling(a/b) par(mfrow=c(b,g)) num<-rep(0,a) sub.plotter<-function(i,fram,framvec,form1,form2,form3,type, xlm,ylm,var1,var2){ temp.i<-fram[framvec <=(i*.10),] plot(form1, data=temp.i, xlim=c(0,xlm), ylim=c(0,ylm), main=((i-1)*.10)) if(type==1){ mod<-lm(form2,data=temp.i) r2<-summary(mod)$adj.r.squared num<-r2 legend("bottomright", legend=signif(r2), col="black") abline(mod) num} else{ if(type==2){ try(mod<-nls(form3, data=temp.i, start=list(x=sx,n=sn), na.action="na.omit"), silent=TRUE) try(x1<-summary(mod)$coefficients[1,1], silent=TRUE) try(n1<-summary(mod)$coefficients[2,1], silent=TRUE) try(lines((1/exp(c(0:70)*x1)*n1)), silent=TRUE) try(num[i]<-AIC(mod), silent=TRUE) try(legend("bottomright", legend=round(num[i],3) , col="black"), silent=TRUE) try((num), silent=TRUE) } }} for(i in 0:a+1){ num<-sub.plotter(i,fram,framvec,form1,form2,form3,type,xlm,ylm) } plot.cor<-function(x){ temp<-a+1 lengthx<-c(1:temp) plot(x~c(1:temp)) m2<-lm(x~c(1:temp)) abline(m2) n<-summary(m2)$adj.r.squared legend("bottomright", legend=signif(n), col="black") slope<-(coef(m2)[2])# slope values<-(num)#values for aic or adj r2 r2ofr2<-(n) #r2 of r2 or AIC output<-data.frame(lengthx,slope,values,r2ofr2) } plot.cor(num) write.csv(plot.cor(num)$output,"output.csv") # can't seem to use paste(substitute(form3),".csv",sep="") to name it at the moment par(mfrow=c(1,1)) } Ben [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.