r - 函数混合模型 B (GAMM)
问题描述
有人用过这个功能吗
“混合.model.B”
? 如果是,我可以找到哪个包?
这是在data.zip部分,数据是白血病.txt 。注意:我正在使用指定的库:
library(nlme)
library(ggplot2)
library(GGally)
library(splines)
library(nlme)
library(fields)
library(lattice)
require(ISLR)
library(grid)
library(dplyr)
library(MASS)
library(mgcv)
library(latticeExtra)
library(fields)
请看下面的代码:
attach(leukemia)
X=model.matrix(height~factor(treatment)*age)
treatment=factor(treatment)
MM=mixed.model.B(age,min(age)-0.5,max(age)+0.5,40,3,2,type="Eilers")
Z=MM[[2]]
Id=factor(rep(1,length(height)))
Z.block4=list(treatment=pdIdent(~Z-1),case=pdSymm(~age))
data.fr <- groupedData(height ~ X[,-1] | Id, data = data.frame(height,X,Z,case,age))
model4 <- lme(height~X[,-1],data=data.fr,random=Z.block4)
## Fitted individual trends for the smooth random intercept and slope model by treatment
b4 <- xyplot(fitted(fit4.gamm$lme) ~ age|factor(treatment),groups=case,col=tim.colors(length(unique(leukemia$case))),
lwd=1,pch=19,data=leukemia,main="Treatment",cex=.35,type="a")
a2 + as.layer(b4)
解决方案
在您列出的链接中,有一个下载 R 文件的链接:http: //halweb.uc3m.es/esp/Personal/personas/durban/esp/web/cursos/Maringa/gam-markdown/Gams-code。 R
在第 349 行,作者构建了函数,我将其复制到下面:
# Mixed model representation of B-splines (using svd)
mixed.model.B<-function(x,xl,xr,ndx,bdeg,pord,type="Eilers"){
Bbasis=bspline(x,xl,xr,ndx,bdeg)
B <- Bbasis$B
m=ncol(B)
D=diff(diag(m),differences=pord)
if(type=="Eilers"){
Z <- B%*%t(D)%*%solve(D%*%t(D))
}else if(type=="SVD"){ print("SVD method")
P.svd=svd(t(D)%*%D)
U=(P.svd$u)[,1:(m-pord)]
d=(P.svd$d)[1:(m-pord)]
Delta=diag(1/sqrt(d))
Z=B%*%U%*%Delta
}
X=NULL
for(i in 0:(pord-1)){
X=cbind(X,x^i)
}
output <- list(X=X,Z=Z)
return(output)
}
推荐阅读
- python - PyQt4 QSqlDatabase:未加载 QMYSQL 驱动程序
- c++ - 在 dlmopen 之后使用 dlsym 时出现分段错误,但使用 dlopen dlsym 工作正常
- react-native - 为什么我使用 `useRef` 创建的我的参考会产生此错误消息?
- google-apps-script - LanguageApp.translate 限制
- demo - “The Mother of All Demos”是用什么语言写的?
- android - Android Studio 警告“此 IDE 正在低内存中运行,而且 android studio 一直没有响应
- json - 为什么 Spring API 在同一个请求中返回 2 个响应?
- python - Python_打印函数的结果
- python - 如何获得 2 个 django 表的每个数量的差异?
- java - 多个循环时的 Java 8 Stream 帮助