遗传力与BLUP.R

library(lme4)

## 载入需要的程辑包：Matrix

data=read.csv("R_ASI遗传力.csv",header = T)
head(data)

##   Line Loc ASI
## 1    1  ZZ   2
## 2    2  ZZ   4
## 3    3  ZZ   4
## 4    4  ZZ   3
## 5    5  ZZ   0
## 6    6  ZZ   5

data$lines=factor(data$Line)
data$env=factor(data$Loc)
data$y=as.numeric(data$ASI)
blp=lmer(y~(1|env)+(1|lines),data=data)
summary(blp)

## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ (1 | env) + (1 | lines)
##    Data: data
## 
## REML criterion at convergence: 2903.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.5108 -0.6603 -0.0469 -0.0081  6.2275 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  lines    (Intercept) 0.03013  0.1736  
##  env      (Intercept) 0.31821  0.5641  
##  Residual             2.48336  1.5759  
## Number of obs: 769, groups:  lines, 203; env, 4
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)    1.877      0.288   6.517

blups= ranef(blp)
names(blups)

## [1] "lines" "env"

lines=blups$lines+blp@beta
res=data.frame(id=rownames(lines),blup=lines)
write.table(res,file="data_blup_result.txt",row.names = F,quote = F,sep="\t")
hist(lines[,1],col="#0AB3CA",border="white",xlab="BLUP of ASI",main="")