It's difficult to see where your error is when you provide no accompanying
data for us to test your code out on (X? Y? Stand?).
However, it looks like you are making your code more complex than it needs
to be. Are you simply trying to a fit a separate linear regression to each
subset of your data identified by a unique standID? If so, you can do so
like this:
# a fake version of your data
mydata <- data.frame(standID=rep(1:5, rep(20, 5)), X=rnorm(100),
Y=rnorm(100))
# split the data up by stand ID
standlist <- split(mydata, mydata$standID)
# fit a regression to each stand and save some results
standEst1 <- lapply(standlist, function(df) {
fit <- summary(lm(Y ~ X, data=df))
co <- coef(fit)
intercept <- co[1, 1]
slope <- co[2, 1]
pValue <- co[2, 4]
mse <- fit$sigma^2
rSquare <- fit$r.squared
c(standID=df$standID[1], intercept=intercept, slope=slope, mse=mse,
rSquare=rSquare, pValue=pValue)
})
# convert the results to a data frame
standEstimates <- data.frame(do.call(rbind, standEst1))
Jean
On Thu, May 16, 2013 at 5:17 AM, rishard <or...@uclive.ac.nz> wrote:
> Hey I'm not really sure what I should put on here, but I am having trouble
> with my R code. I am trying to get the p-values, R^2s etc for a number of
> different groups of variables that are all in one dataset.
>
> This is the code:
>
> #Stand counter
> st<-1
> #Collections
> stands<-numeric(67)
> slopes<-numeric(67)
> intercepts<-numeric(67)
> mses<-numeric(67)
> rsquares<-numeric(67)
> pValues<-numeric(67)
> #Start lists for X and Y values within each stand
> xi<-numeric(0)
> yi<-numeric(0)
> #Set the first element to the starting X and Y values
> xi[1]=X[1]
> yi[1]=Y[1]
> #Start looping working through your data, record by record
> for (i in 2:length(X)) {
> #If you are in the same stand as on the last record, continue to
> #collect X and Y values
> if(Stand[i]==Stand[i-1]) {
> xi=cbind(xi,X[i])
> yi=cbind(yi,Y[i])
> } else {
> #If a new stand is encountered make your linear model and
> #collect statistics
> model<-lm(yi~xi)
> stands[st]<-Stand[i-1]
> intercepts[st]<-model$coefficients[1]
> slopes[st]<-model$coefficients[2]
> mses[st]<-sum(resid(model)^2)/(length(xi)-2)
> ssr<-var(yi)*(length(xi)-1)-sum(resid(model)^2)
> rsquares[st]<-ssr/(var(yi)*(length(xi)-1))
> fRatio<-ssr/mses[st]
> pValues[st]<-1-pf(fRatio,1,length(xi)-2)
> #Increment the stand number, zero the within stand collections,
> #and start again
> st<-st+1
> xi<-numeric(0)
> yi<-numeric(0)
> xi[1]=X[i]
> yi[1]=Y[i]
> }
> }
> #Make your data set
>
> standEstimates<-data.frame(standID=stands,intercept=intercepts,slop=slopes,mse=mses,rSquare=rsquares,pValue=pValues)
>
> The standEstimate outputs look like this:
>
> standID intercept slop mse rSquare pValue
> 1 6833 319.2470 NA 0 NA NA
> 2 756 708.7470 NA 0 NA NA
> 3 795 508.2290 NA 0 NA NA
> 4 1249 503.1460 NA 0 NA NA
> 5 1331 703.0620 NA 0 NA NA
> 6 1417 747.7620 NA 0 NA NA
> 7 4715 549.3400 NA 0 NA NA
> 8 4850 603.9940 NA 0 NA NA
> 9 2105 573.3270 NA 0 NA NA
> Etc etc
>
> and I get these warnings:
>
> 1: In rsquares[st] <- ssr/(var(yi) * (length(xi) - 1)) :
> number of items to replace is not a multiple of replacement length
> 2: In pValues[st] <- 1 - pf(fRatio, 1, length(xi) - 2) :
> number of items to replace is not a multiple of replacement length
> 3: In rsquares[st] <- ssr/(var(yi) * (length(xi) - 1)) :
> number of items to replace is not a multiple of replacement length
> 4: In pValues[st] <- 1 - pf(fRatio, 1, length(xi) - 2) :
> number of items to replace is not a multiple of replacement length
> 5: In rsquares[st] <- ssr/(var(yi) * (length(xi) - 1)) :
> number of items to replace is not a multiple of replacement length
> 6: In pValues[st] <- 1 - pf(fRatio, 1, length(xi) - 2) :
> number of items to replace is not a multiple of replacement length
> 7: In rsquares[st] <- ssr/(var(yi) * (length(xi) - 1)) :
> number of items to replace is not a multiple of replacement length
> 8: In pValues[st] <- 1 - pf(fRatio, 1, length(xi) - 2) :
> number of items to replace is not a multiple of replacement length
> 9: In rsquares[st] <- ssr/(var(yi) * (length(xi) - 1)) :
> number of items to replace is not a multiple of replacement length
> 10: In pValues[st] <- 1 - pf(fRatio, 1, length(xi) - 2) :
> number of items to replace is not a multiple of replacement length
> 11: In rsquares[st] <- ssr/(var(yi) * (length(xi) - 1)) :
> number of items to replace is not a multiple of replacement length
> 12: In pValues[st] <- 1 - pf(fRatio, 1, length(xi) - 2) :
> number of items to replace is not a multiple of replacement length
> 13: In rsquares[st] <- ssr/(var(yi) * (length(xi) - 1)) :
> number of items to replace is not a multiple of replacement length
> 14: In pValues[st] <- 1 - pf(fRatio, 1, length(xi) - 2) :
> etc etc
>
> Sorry if I haven't set this post out right, or haven't provided enough
> information. But can anyone see why it is not giving me any returns for
> R^2
> etc?
>
>
>
> --
> View this message in context:
> http://r.789695.n4.nabble.com/R-looping-help-tp4667180.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> 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.
>
[[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.
