residual of LM

怎样做多度与环境的线性回归好呢？(未解决)

Ntotal = 50 *50
xi.1 <- runif(n=Ntotal, min = 0,max = 500)
yi.1 <- runif(n=Ntotal, min = 0,max = 300)
plot(xi.1, yi.1)

#
Nparents= 50 ; Nkids=50; SD=3
x0 <- runif(n=Nparents, min = 10, max = 490)
y0 <- runif(n=Nparents, min = 10, max = 290)
xi.c <- rep(x0, each=Nkids) + rnorm(n=Nparents * Nkids, mean=0, sd=SD)
yi.c <- rep(y0, each=Nkids) + rnorm(n=Nparents * Nkids, mean=0, sd=SD)
plot(xi.c, yi.c)

#
xf <- factor(ceiling(xi.c/10), levels=1:50)
yf <- factor(ceiling(yi.c/10), levels=1:30)
xyTab.c <- table(yf, xf) 
image(t(xyTab.c))
points(xi.c/500, yi.c/300, pch='*')

#
xi <- c(xi.1, xi.c); yi <- c(yi.1, yi.c)
plot(xi, yi)

#
xif <- factor(ceiling(xi/10),levels=1:50)
yif <- factor(ceiling(yi/10),levels=1:30)
xyi.tab <- table(yif, xif)
#
y <- c(xyi.tab)
x <- c(xyTab.c)
#
plot(x,y)

plot(x,log(y))

plot(log(x),log(y))

#
xlog <- log1p(x)
ylog <- log1p(y)
#
lm0 <- lm(y~x)
lm1 <- lm(ylog~xlog)
#
AIC(lm0)

## [1] 4984.153

AIC(lm1)

## [1] 2060.851

summary(lm0)

## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.1246 -0.6473 -0.6473  0.3527  5.3527 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.647263   0.034149   48.24   <2e-16 ***
## x           1.011642   0.005577  181.41   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.273 on 1498 degrees of freedom
## Multiple R-squared:  0.9565, Adjusted R-squared:  0.9564 
## F-statistic: 3.291e+04 on 1 and 1498 DF,  p-value: < 2.2e-16

summary(lm1)

## 
## Call:
## lm(formula = ylog ~ xlog)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.84862 -0.15547 -0.04764  0.24999  1.23082 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.84862    0.01322   64.19   <2e-16 ***
## xlog         0.72679    0.01507   48.23   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4803 on 1498 degrees of freedom
## Multiple R-squared:  0.6083, Adjusted R-squared:  0.608 
## F-statistic:  2326 on 1 and 1498 DF,  p-value: < 2.2e-16

#
y0.res <- resid(lm0)
y0.fit <- predict(lm0)
plot(y0.fit, y0.res)

plot(lm0)

#
y1.res <- resid(lm1)
y1.fit <- predict(lm1)
plot(y1.fit, y1.res)

plot(lm1)

#
plot(lm0,2)

plot(lm1,2)

#
hist(log(y))

#
fit0 <- lm(exp(y) ~ x+I(x^2))
plot(fit0)

#
lambda0 <-  predict(fit0) 
#
plot(c(xyTab.c), c(xyi.tab))

plot(fit0)