R: Converting do.call()-summary in summary -


as stepaic() function mass package has problems when used within function, use do.call() (described here). problem sounds easy could't find solution it: when use do.call() lm() model several raster layers, layers saved within model. if want print summary() of model, writes layers in output , gets confusing. how "normal" summary output, without using do.call?

here short example:

create list of raster layers:

xz.list  <- lapply(1:5,function(x){   r1 <- raster(ncol=3, nrow=3)   values(r1) <- 1:ncell(r1)   r1 }) 

convert them in data.frame:

xz<-getvalues(stack(xz.list))  xz <- as.data.frame(xz) 

use do.call lm model:

fit1<-do.call("lm", list(xz[,1] ~ . , data = xz)) 

the summary() output looks this:

summary(fit1)  call: lm(formula = xz[, 1] ~ ., data = structure(list(layer.1 = 1:9,      layer.2 = 1:9, layer.3 = 1:9, layer.4 = 1:9, layer.5 = 1:9), .names = c("layer.1",  "layer.2", "layer.3", "layer.4", "layer.5"), row.names = c(na,  -9l), class = "data.frame"))  residuals:        min         1q     median         3q        max  -9.006e-16 -2.472e-16 -2.031e-16 -1.370e-16  1.724e-15   coefficients: (4 not defined because of singularities)              estimate std. error   t value pr(>|t|)     (intercept) 1.184e-15  5.784e-16 2.047e+00   0.0798 .   layer.1     1.000e+00  1.028e-16 9.729e+15   <2e-16 *** layer.2            na         na        na       na     layer.3            na         na        na       na     layer.4            na         na        na       na     layer.5            na         na        na       na     --- signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1   residual standard error: 7.962e-16 on 7 degrees of freedom multiple r-squared:     1,  adjusted r-squared:     1  f-statistic: 9.465e+31 on 1 , 7 df,  p-value: < 2.2e-16  

this doesn't bad in small example, becomes mess when using 10 or more raster layers 32k values each. make output use summary(lm) function without do.call:

fit<-lm(xz[,1] ~ . , data=xz) summary(fit)  call: lm(formula = xz[, 1] ~ ., data = xz)  residuals:        min         1q     median         3q        max  -9.006e-16 -2.472e-16 -2.031e-16 -1.370e-16  1.724e-15   coefficients: (4 not defined because of singularities)              estimate std. error   t value pr(>|t|)     (intercept) 1.184e-15  5.784e-16 2.047e+00   0.0798 .   layer.1     1.000e+00  1.028e-16 9.729e+15   <2e-16 *** layer.2            na         na        na       na     layer.3            na         na        na       na     layer.4            na         na        na       na     layer.5            na         na        na       na     --- signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1   residual standard error: 7.962e-16 on 7 degrees of freedom multiple r-squared:     1,  adjusted r-squared:     1  f-statistic: 9.465e+31 on 1 , 7 df,  p-value: < 2.2e-16  

you can redefine lm function :

lm <- function(form, ...) { fm <- stats::lm(form,...);                              fm$call <- form; fm } 

testing :

fit2<-do.call("lm", list(xz[,1] ~ . , data = xz))  summary(fit2)  call: xz[, 1] ~ .  residuals:        min         1q     median         3q        max  -9.006e-16 -2.472e-16 -2.031e-16 -1.370e-16  1.724e-15   coefficients: (4 not defined because of singularities)              estimate std. error   t value pr(>|t|)     (intercept) 1.184e-15  5.784e-16 2.047e+00   0.0798 .   layer.1     1.000e+00  1.028e-16 9.729e+15   <2e-16 *** layer.2            na         na        na       na     layer.3            na         na        na       na     layer.4            na         na        na       na     layer.5            na         na        na       na     --- signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1  residual standard error: 7.962e-16 on 7 degrees of freedom multiple r-squared:      1, adjusted r-squared:      1  f-statistic: 9.465e+31 on 1 , 7 df,  p-value: < 2.2e-16 

Comments

Popular posts from this blog

javascript - Count length of each class -

What design pattern is this code in Javascript? -

hadoop - Restrict secondarynamenode to be installed and run on any other node in the cluster -