h2_mzdz | R Documentation |
Heritability estimation according to twin correlations
h2_mzdz(
mzDat = NULL,
dzDat = NULL,
rmz = NULL,
rdz = NULL,
nmz = NULL,
ndz = NULL,
selV = NULL
)
mzDat |
a data frame for monzygotic twins (MZ). |
dzDat |
a data frame for dizygotic twins (DZ). |
rmz |
correlation for MZ twins. |
rdz |
correlation for DZ twins. |
nmz |
sample size for MZ twins. |
ndz |
sample size for DZ twins. |
selV |
names of variables for twin and cotwin. |
Given MZ/DZ data or their correlations and sample sizes, it obtains heritability and variance estimates under an ACE model as in \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1038/s41562-023-01530-y")} and \insertCitekeeping95;textualgap.
A data.frame with variables h2, c2, e2, vh2, vc2, ve2.
elks12gap
\insertAllCited## Not run:
library(mvtnorm)
set.seed(12345)
mzm <- as.data.frame(rmvnorm(195, c(22.75,22.75),
matrix(2.66^2*c(1, 0.67, 0.67, 1), 2)))
dzm <- as.data.frame(rmvnorm(130, c(23.44,23.44),
matrix(2.75^2*c(1, 0.32, 0.32, 1), 2)))
mzw <- as.data.frame(rmvnorm(384, c(21.44,21.44),
matrix(3.08^2*c(1, 0.72, 0.72, 1), 2)))
dzw <- as.data.frame(rmvnorm(243, c(21.72,21.72),
matrix(3.12^2*c(1, 0.33, 0.33, 1), 2)))
selVars <- c('bmi1','bmi2')
names(mzm) <- names(dzm) <- names(mzw) <- names(dzw) <- selVars
ACE_CI <- function(mzData,dzData,n.sim=5,selV=NULL,verbose=TRUE)
{
ACE_obs <- h2_mzdz(mzDat=mzData,dzDat=dzData,selV=selV)
cat("\n\nheritability according to correlations\n\n")
print(format(ACE_obs,digits=3),row.names=FALSE)
nmz <- nrow(mzData)
ndz <- nrow(dzData)
r <- data.frame()
for(i in 1:n.sim)
{
cat("\rRunning # ",i,"/", n.sim,"\r",sep="")
sampled_mz <- sample(1:nmz, replace=TRUE)
sampled_dz <- sample(1:ndz, replace=TRUE)
mzDat <- mzData[sampled_mz,]
dzDat <- dzData[sampled_dz,]
ACE_i <- h2_mzdz(mzDat=mzDat,dzDat=dzDat,selV=selV)
if (verbose) print(ACE_i)
r <- rbind(r,ACE_i)
}
m <- apply(r,2,mean,na.rm=TRUE)
s <- apply(r,2,sd,na.rm=TRUE)
allr <- data.frame(mean=m,sd=s,lcl=m-1.96*s,ucl=m+1.96*s)
print(format(allr,digits=3))
}
ACE_CI(mzm,dzm,n.sim=500,selV=selVars,verbose=FALSE)
ACE_CI(mzw,dzw,n.sim=500,selV=selVars,verbose=FALSE)
## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.