Analyzing twin data with 'mets'

Table of Contents

Installation

Install dependencies (R>=2.15) :

install.packages(c("mets","cmprsk"), dependencies=TRUE)

OBS: At this point you might have to restart R to flush the cache of previously installed versions of the packages. If you have previously installed timereg and lava, make sure that you have the current versions installed (timereg: 1.8.4, lava: 1.2.6).

Load simulated data

library(mets)

The dataset prt contains (simulated) observations on prostate cancer with the following columns

country
Country (Denmark,Finland,Norway,Sweden)
time
exit time (censoring,death or prostate cancer)
status
Status (censoring=0,death=1 or prostate cancer=2)
zyg
Zygosity (DZ,MZ)
id
Twin id number
cancer
cancer indicator (status=2)
data(prt)
head(prt)
   country     time status zyg id cancer
31 Denmark 96.98833      1  DZ  1      0
32 Denmark 80.88885      1  DZ  1      0
39 Denmark 68.04498      1  DZ  3      0
40 Denmark 61.45903      1  DZ  3      0
51 Denmark 78.78068      1  DZ  5      0
52 Denmark 90.36252      1  DZ  5      0

Status table

prtwide <- fast.reshape(prt,id="id")
ftable(status1~status2,prtwide)
        status1    0    1    2
status2                       
0               9278  883  156
1                936 2308  193
2                163  199  106

Estimation of cumulative incidence

times <- seq(40,100,by=2)
cifmod <- comp.risk(Hist(time,status)~+1+cluster(id),data=prt,
                    cause=2,n.sim=0,
                    times=times,conservative=1,max.clust=NULL,model="fg")

theta.des <- model.matrix(~-1+factor(zyg),data=prt) ## design for MZ/DZ status
or1 <- or.cif(cifmod,data=prt,cause1=2,cause2=2,theta.des=theta.des,
              score.method="fisher.scoring",same.cens=TRUE)
summary(or1)
or1$score
OR for dependence for competing risks

OR of cumulative incidence for cause1= 2  and cause2= 2
              log-ratio Coef.    SE    z    P-val Ratio    SE
factor(zyg)DZ           0.785 0.221 3.55 3.82e-04  2.19 0.485
factor(zyg)MZ           2.100 0.278 7.56 4.11e-14  8.14 2.260
              [,1]
[1,] -5.179525e-09
[2,]  4.202363e-08
pcif <- predict(cifmod,X=1,resample.iid=0,uniform=0,se=0)
plot(pcif,multiple=1,se=0,uniform=0,ylim=c(0,0.15))

pcif.png

Assumes that the censoring of the two twins are independent (when they are the same):

incorrect.or1 <- or.cif(cifmod,data=prt,cause1=2,cause2=2,theta.des=theta.des, 
                        theta=c(2.8,8.6),score.method="fisher.scoring")
summary(incorrect.or1)
## not  bad
incorrect.or1$score
OR for dependence for competing risks

OR of cumulative incidence for cause1= 2  and cause2= 2
              log-ratio Coef.    SE    z    P-val  Ratio       SE
factor(zyg)DZ            2.84 0.469 6.06 1.37e-09   17.1     8.04
factor(zyg)MZ            8.61 5.810 1.48 1.39e-01 5460.0 31800.00

Correcting for country

table(prt$country)

times <- seq(40,100,by=2)
cifmodl <-comp.risk(Hist(time,status)~-1+factor(country)+cluster(id),data=prt,
                    cause=2,n.sim=0,times=times,conservative=1,
                    max.clust=NULL,cens.model="aalen")
pcifl <- predict(cifmodl,X=diag(4),se=0,uniform=0)
plot(pcifl,multiple=1,se=0,uniform=0,col=1:4,ylim=c(0,0.2))
legend("topleft",levels(prt$country),col=1:4,lty=1)

pcifl.png

Design for MZ/DZ status

theta.des <- model.matrix(~-1+factor(zyg),data=prt) 
or.country <- or.cif(cifmodl,data=prt,cause1=2,cause2=2,theta.des=theta.des,
                     theta=c(0.8,2.1),score.method="fisher.scoring",same.cens=TRUE)

summary(or.country)
OR for dependence for competing risks

OR of cumulative incidence for cause1= 2  and cause2= 2
              log-ratio Coef.    SE    z    P-val Ratio    SE
factor(zyg)DZ           0.736 0.234 3.15 1.66e-03  2.09 0.488
factor(zyg)MZ           1.860 0.279 6.67 2.54e-11  6.44 1.800

Concordance estimation

Ignoring country. Computing casewise, using prodlim. CIF:

outm <- prodlim(Hist(time,status)~+1,data=prt)

times <- 60:100
## cause is 2 (second cause)
cifmz <- predict(outm,cause=2,time=times,newdata=data.frame(zyg="MZ"))
cifdz <- predict(outm,cause=2,time=times,newdata=data.frame(zyg="DZ"))
### casewise 
pp33 <- bicomprisk(Hist(time,status)~strata(zyg)+id(id),data=prt,cause=c(2,2),prodlim=TRUE)
pp33dz <- pp33$model$"DZ"
pp33mz <- pp33$model$"MZ"
concdz <- predict(pp33dz,cause=1,time=times,newdata=data.frame(zyg="DZ"))
concmz <- predict(pp33mz,cause=1,time=times,newdata=data.frame(zyg="MZ"))
par(mfrow=c(1,2))
plot(times,concdz,ylim=c(0,0.1),type="s")
lines(pcif$time,pcif$P1^2,col=2)
title(main="DZ Conc. Prostate cancer")
plot(times,concmz,ylim=c(0,0.1),type="s")
title(main="MZ Conc. Prostate cancer")
lines(pcif$time,pcif$P1^2,col=2)

concordance.png

par(mfrow=c(1,1))
cdz <- casewise(pp33dz,outm,cause.marg=2)
cmz <- casewise(pp33mz,outm,cause.marg=2)             
plot(cmz,ci=NULL,ylim=c(0,0.5),xlim=c(60,100),legend=TRUE,col=c(3,2,1))
par(new=TRUE)
plot(cdz,ci=NULL,ylim=c(0,0.5),xlim=c(60,100),legend=TRUE)

casewisea.png

Similar analyses using comp.risk for competing risks leads to tests for equal concordance and more correct standard errors

p33 <- bicomprisk(Hist(time,status)~strata(zyg)+id(id),data=prt,cause=c(2,2),return.data=1)

p33dz <- p33$model$"DZ"$comp.risk
p33mz <- p33$model$"MZ"$comp.risk
head(cbind(p33mz$time, p33mz$P1, p33mz$se.P1))
head(cbind(p33dz$time, p33dz$P1, p33dz$se.P1))
         [,1]        [,2]         [,3]
[1,] 60.88384 0.001354486 0.0006759148
[2,] 64.98252 0.001738665 0.0007767791
[3,] 66.34227 0.002145175 0.0008759241
[4,] 67.23626 0.002553690 0.0009656368
[5,] 67.96152 0.002980112 0.0010544136
[6,] 68.37310 0.003852670 0.0012192761
         [,1]         [,2]         [,3]
[1,] 58.85519 0.0001741916 0.0001740997
[2,] 67.87387 0.0004044091 0.0002883926
[3,] 69.55123 0.0006488647 0.0003777479
[4,] 70.83183 0.0009069944 0.0004570724
[5,] 71.05738 0.0011672691 0.0005255212
[6,] 71.06602 0.0014276382 0.0005859026

Test for genetic effect, needs other form of bicomprisk with iid decomp

conc1 <- p33dz
conc2 <- p33mz

test.conc(p33dz,p33mz);
$test
           cum dif.         sd        z        pval
pepe-mori 0.3937372 0.09841628 4.000732 6.31468e-05

$mintime
[1] 60.88384

$maxtime
[1] 96.92463

$same.cluster
[1] FALSE

attr(,"class")
[1] "testconc"

OR expression of difference in concordance functions and Gray test

data33mz <- p33$model$"MZ"$data
data33mz$zyg <- 1
data33dz <- p33$model$"DZ"$data
data33dz$zyg <- 0
data33 <- rbind(data33mz,data33dz)

library(cmprsk)
ftime <- data33$time
fstatus <- data33$status
table(fstatus)
fstatus
   0    1    2 
9597  106 4519
group <- data33$zyg
graytest <- cuminc(ftime,fstatus,group)
graytest
Tests:
      stat           pv df
1 28.82416 7.925617e-08  1
2 33.79236 6.131919e-09  1
Estimates and Variances:
$est
              20         40           60          80        100
0 1 0.0000000000 0.00000000 0.0001741916 0.006741025 0.01880244
1 1 0.0000000000 0.00000000 0.0006710172 0.017420360 0.05031415
0 2 0.0006970762 0.01974882 0.1141800067 0.504364854 0.93797293
1 2 0.0009363302 0.01655314 0.0948098327 0.443996722 0.90692430

$var
              20           40           60           80          100
0 1 0.000000e+00 0.000000e+00 3.034323e-08 2.115863e-06 9.493584e-06
1 1 0.000000e+00 0.000000e+00 2.250627e-07 9.173278e-06 5.102841e-05
0 2 8.094463e-08 2.487399e-06 1.556735e-05 6.990685e-05 4.769058e-05
1 2 1.752378e-07 3.424511e-06 2.388136e-05 1.271394e-04 1.171775e-04
zygeffect <- comp.risk(Hist(time,status)~const(zyg),
                  data=data33,cause=1,
                  cens.model="aalen",model="logistic",conservative=1)
summary(zygeffect)
Competing risks Model 

Test for nonparametric terms 

Test for non-significant effects 
            Supremum-test of significance p-value H_0: B(t)=0
(Intercept)                          26.8                   0

Test for time invariant effects 
                  Kolmogorov-Smirnov test p-value H_0:constant effect
(Intercept)                          2.22                           0
                    Cramer von Mises test p-value H_0:constant effect
(Intercept)                          36.3                           0

Parametric terms : 
           Coef.    SE Robust SE    z    P-val
const(zyg) 0.944 0.218     0.218 4.34 1.45e-05
   
  Call: 
comp.risk(Hist(time, status) ~ const(zyg), data = data33, cause = 1, 
    cens.model = "aalen", model = "logistic", conservative = 1)

Liability model, ignoring censoring

(M <- with(prt, table(cancer,zyg)))
      zyg
cancer    DZ    MZ
     0 17408 10872
     1   583   359
coef(lm(cancer~-1+zyg,prt))
     zygDZ      zygMZ 
0.03240509 0.03196510

Saturated model

bpmz <- biprobit(cancer~1 + cluster(id), 
             data=subset(prt,zyg=="MZ"), eqmarg=TRUE)

logLik(bpmz) # Log-likelihood
AIC(bpmz) # AIC
coef(bpmz) # Parameter estimates
vcov(bpmz) # Asymptotic covariance
summary(bpmz) # concordance, case-wise, tetrachoric correlations, ...
'log Lik.' -1472.972 (df=2)
[1] 2949.943
(Intercept)  atanh(rho) 
 -1.8539454   0.8756506
             (Intercept)   atanh(rho)
(Intercept) 0.0007089726 0.0003033296
atanh(rho)  0.0003033296 0.0044023587

              Estimate    Std.Err          Z p-value
(Intercept)  -1.853945   0.026627 -69.627727       0
atanh(rho)    0.875651   0.066350  13.197393       0

    n pairs 
11231  5473 
Score: -3.453e-05 5.123e-06
logLik: -1472.972 
Variance of latent residual term = 1 (standard probit link) 

                        Estimate 2.5%    97.5%  
Tetrachoric correlation 0.70423  0.63252 0.76398
Concordance             0.01131  0.00886 0.01443
Case-wise/Conditional   0.35487  0.29391 0.42094
Marginal                0.03187  0.02834 0.03583
bp0 <- biprobit(cancer~1 + cluster(id)+strata(zyg), data=prt)
summary(bp0)
------------------------------------------------------------
Strata 'DZ'

              Estimate    Std.Err          Z p-value
(Intercept)  -1.846841   0.019247 -95.955243       0
atanh(rho)    0.418065   0.050421   8.291446       0

    n pairs 
17991  8749 
Score: -0.001842 -0.0006881
logLik: -2536.242 
Variance of latent residual term = 1 (standard probit link) 

                        Estimate 2.5%    97.5%  
Tetrachoric correlation 0.39530  0.30882 0.47529
Concordance             0.00486  0.00361 0.00655
Case-wise/Conditional   0.15019  0.11459 0.19443
Marginal                0.03239  0.02976 0.03523

------------------------------------------------------------
Strata 'MZ'

              Estimate    Std.Err          Z p-value
(Intercept)  -1.853945   0.026627 -69.627727       0
atanh(rho)    0.875651   0.066350  13.197393       0

    n pairs 
11231  5473 
Score: -3.453e-05 5.123e-06
logLik: -1472.972 
Variance of latent residual term = 1 (standard probit link) 

                        Estimate 2.5%    97.5%  
Tetrachoric correlation 0.70423  0.63252 0.76398
Concordance             0.01131  0.00886 0.01443
Case-wise/Conditional   0.35487  0.29391 0.42094
Marginal                0.03187  0.02834 0.03583

Equal marginals MZ/DZ

bp1 <- bptwin(cancer~1,zyg="zyg",DZ="DZ",id="id",type="u",data=prt)
(s <- summary(bp1))
                 Estimate     Std.Err           Z p-value
(Intercept)     -1.849284    0.015601 -118.539777       0
atanh(rho) MZ    0.877667    0.065815   13.335456       0
atanh(rho) DZ    0.417475    0.050276    8.303615       0

 MZ/DZ       Complete pairs MZ/DZ
 11231/17991 5473/8749           

                           Estimate 2.5%    97.5%  
Tetrachoric correlation MZ 0.70525  0.63436 0.76438
Tetrachoric correlation DZ 0.39480  0.30854 0.47462

MZ:
                     Estimate 2.5%     97.5%   
Concordance           0.01149  0.00942  0.01400
Casewise Concordance  0.35672  0.29764  0.42049
Marginal              0.03221  0.03007  0.03449
Rel.Recur.Risk       11.07524  9.15861 12.99187
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.00482  0.00363 0.00640
Casewise Concordance 0.14956  0.11441 0.19315
Marginal             0.03221  0.03007 0.03449
Rel.Recur.Risk       4.64343  3.44806 5.83880

                         Estimate 2.5%    97.5%  
Broad-sense heritability 0.62090  0.41075 0.83104

Components (concordance,cor,…) can be extracted from returned list

s$all
                               Estimate        2.5%        97.5%
Broad-sense heritability    0.620895137 0.410750804  0.831039470
Tetrachoric correlation MZ  0.705248651 0.634356556  0.764377527
Tetrachoric correlation DZ  0.394801083 0.308543835  0.474618270
MZ Concordance              0.011489242 0.009421632  0.014004180
MZ Casewise Concordance     0.356715720 0.297643978  0.420492296
MZ Marginal                 0.032208397 0.030073567  0.034489384
MZ Rel.Recur.Risk          11.075239652 9.158610651 12.991868652
DZ Concordance              0.004817009 0.003625030  0.006398416
DZ Casewise Concordance     0.149557550 0.114405842  0.193154111
DZ Marginal                 0.032208397 0.030073567  0.034489384
DZ Rel.Recur.Risk           4.643433455 3.448063061  5.838803849

Likelihood Ratio Test

compare(bp0,bp1)
	- Likelihood ratio test -

data:  
chisq = 0.0468, df = 1, p-value = 0.8288
sample estimates:
log likelihood (model 1) log likelihood (model 2) 
               -4009.213                -4009.237

Polygenic Libability model via te bptwin function (type can be a subset of "acde", or "flex" for stratitified, "u" for random effects model with same marginals for MZ and DZ)

bp2 <- bptwin(cancer~1,zyg="zyg",DZ="DZ",id="id",type="ace",data=prt)
summary(bp2)
             Estimate   Std.Err         Z p-value
(Intercept)  -3.40624   0.19032 -17.89736  0.0000
log(var(A))   0.74503   0.25710   2.89787  0.0038
log(var(C))  -1.25112   1.04238  -1.20024  0.2300

 MZ/DZ       Complete pairs MZ/DZ
 11231/17991 5473/8749           

                   Estimate 2.5%     97.5%   
A                   0.62090  0.41075  0.83104
C                   0.08435 -0.09373  0.26244
E                   0.29475  0.22992  0.35959
MZ Tetrachoric Cor  0.70525  0.63436  0.76438
DZ Tetrachoric Cor  0.39480  0.30854  0.47462

MZ:
                     Estimate 2.5%     97.5%   
Concordance           0.01149  0.00942  0.01400
Casewise Concordance  0.35672  0.29764  0.42049
Marginal              0.03221  0.03007  0.03449
Rel.Recur.Risk       11.07524  9.15861 12.99187
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.00482  0.00363 0.00640
Casewise Concordance 0.14956  0.11441 0.19315
Marginal             0.03221  0.03007 0.03449
Rel.Recur.Risk       4.64343  3.44806 5.83880

                         Estimate 2.5%    97.5%  
Broad-sense heritability 0.62090  0.41075 0.83104

Liability model, Inverse Probability Weighting

Probability weights based on Aalen's additive model

prtw <- ipw(Surv(time,status==0)~country, data=prt,
            cluster="id",weightname="w") 
plot(0,type="n",xlim=range(prtw$time),ylim=c(0,1),xlab="Age",ylab="Probability")
count <- 0
for (l in unique(prtw$country)) {
    count <- count+1
    prtw <- prtw[order(prtw$time),]
    with(subset(prtw,country==l), 
         lines(time,w,col=count,lwd=2))
}
legend("topright",legend=unique(prtw$country),col=1:4,pch=-1,lty=1)

ipw.png

bpmzIPW <- biprobit(cancer~1 + cluster(id), 
                    data=subset(prtw,zyg=="MZ"), 
                    weight="w")
(smz <- summary(bpmzIPW))
              Estimate    Std.Err          Z p-value
(Intercept)  -1.226276   0.043074 -28.469378       0
atanh(rho)    0.912670   0.100316   9.097911       0

    n pairs 
 2722   997 
Score: 3.329e-05 -2.252e-05
logLik: -6703.246 
Variance of latent residual term = 1 (standard probit link) 

                        Estimate 2.5%    97.5%  
Tetrachoric correlation 0.72241  0.61446 0.80381
Concordance             0.05490  0.04221 0.07113
Case-wise/Conditional   0.49887  0.41321 0.58460
Marginal                0.11005  0.09514 0.12696

Comparison with CIF

plot(pcif,multiple=1,se=1,uniform=0,ylim=c(0,0.15))
abline(h=smz$prob["Marginal",],lwd=c(2,1,1))
## Wrong estimates:
abline(h=summary(bpmz)$prob["Marginal",],lwd=c(2,1,1),col="lightgray")

cifMZ.png

Concordance estimates

plot(pp33mz,ylim=c(0,0.1))
abline(h=smz$prob["Concordance",],lwd=c(2,1,1))
## Wrong estimates:
abline(h=summary(bpmz)$prob["Concordance",],lwd=c(2,1,1),col="lightgray")

conc2.png

ACE model with IPW

bp3 <- bptwin(cancer~1,zyg="zyg",DZ="DZ",id="id",
              type="ace",data=prtw,weight="w")
summary(bp3)
             Estimate   Std.Err         Z p-value
(Intercept)  -2.31618   0.18673 -12.40359   0e+00
log(var(A))   0.85390   0.22689   3.76347   2e-04
log(var(C)) -29.43218   4.34216  -6.77823   0e+00

 MZ/DZ     Complete pairs MZ/DZ
 4716/8835 997/1809            

                   Estimate 2.5%    97.5%  
A                  0.70138  0.60824 0.79452
C                  0.00000      NaN     NaN
E                  0.29862  0.20548 0.39176
MZ Tetrachoric Cor 0.70138  0.59586 0.78310
DZ Tetrachoric Cor 0.35069  0.30328 0.39637

MZ:
                     Estimate 2.5%    97.5%  
Concordance          0.04857  0.03963 0.05940
Casewise Concordance 0.47238  0.39356 0.55260
Marginal             0.10281  0.09463 0.11161
Rel.Recur.Risk       4.59457  3.79490 5.39425
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.02515  0.02131 0.02965
Casewise Concordance 0.24461  0.21892 0.27226
Marginal             0.10281  0.09463 0.11161
Rel.Recur.Risk       2.37919  2.13966 2.61872

                         Estimate 2.5%    97.5%  
Broad-sense heritability 0.70138  0.60824 0.79452

Equal marginals but free variance structure between MZ and DZ

bp4 <- bptwin(cancer~1,zyg="zyg",DZ="DZ",id="id",
              type="u",data=prtw,weight="w")
summary(bp4)
                Estimate    Std.Err          Z p-value
(Intercept)    -1.266427   0.024091 -52.568381       0
atanh(rho) MZ   0.898548   0.098841   9.090866       0
atanh(rho) DZ   0.312574   0.073668   4.243006       0

 MZ/DZ     Complete pairs MZ/DZ
 4716/8835 997/1809            

                           Estimate 2.5%    97.5%  
Tetrachoric correlation MZ 0.71559  0.60742 0.79771
Tetrachoric correlation DZ 0.30278  0.16662 0.42760

MZ:
                     Estimate 2.5%    97.5%  
Concordance          0.04974  0.04044 0.06104
Casewise Concordance 0.48442  0.40185 0.56785
Marginal             0.10268  0.09453 0.11144
Rel.Recur.Risk       4.71777  3.88751 5.54802
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.02269  0.01667 0.03081
Casewise Concordance 0.22097  0.16448 0.29013
Marginal             0.10268  0.09453 0.11144
Rel.Recur.Risk       2.15203  1.53917 2.76490

                         Estimate 2.5%    97.5%  
Broad-sense heritability 0.82563  0.50195 1.14931

Check convergence

mean(score(bp4)^2)
[1] 2.723881e-14

Liability model, adjusting for covariates

Main effect of country

bp6 <- bptwin(cancer~country,zyg="zyg",DZ="DZ",id="id",
              type="ace",data=prtw,weight="w")
summary(bp6)
                Estimate   Std.Err         Z p-value
(Intercept)     -2.81553   0.23889 -11.78590   0e+00
countryFinland   0.87558   0.16123   5.43061   0e+00
countryNorway    0.68483   0.17762   3.85567   1e-04
countrySweden    0.77248   0.12350   6.25468   0e+00
log(var(A))      0.77724   0.23186   3.35220   8e-04
log(var(C))    -28.96268   4.10060  -7.06303   0e+00

 MZ/DZ     Complete pairs MZ/DZ
 4716/8835 997/1809            

                   Estimate 2.5%    97.5%  
A                  0.68509  0.58704 0.78313
C                  0.00000      NaN     NaN
E                  0.31491  0.21687 0.41296
MZ Tetrachoric Cor 0.68509  0.57428 0.77124
DZ Tetrachoric Cor 0.34254  0.29262 0.39060

MZ:
                     Estimate 2.5%    97.5%  
Concordance          0.02236  0.01588 0.03141
Casewise Concordance 0.39194  0.30778 0.48305
Marginal             0.05705  0.04654 0.06977
Rel.Recur.Risk       6.86967  5.08343 8.65591
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.00989  0.00700 0.01394
Casewise Concordance 0.17329  0.14505 0.20570
Marginal             0.05705  0.04654 0.06977
Rel.Recur.Risk       3.03735  2.56114 3.51356

                         Estimate 2.5%    97.5%  
Broad-sense heritability 0.68509  0.58704 0.78313

bp7 <- bptwin(cancer~country,zyg="zyg",DZ="DZ",id="id",
              type="u",data=prtw,weight="w")
summary(bp7)
                 Estimate    Std.Err          Z p-value
(Intercept)     -1.581478   0.051318 -30.817030   0e+00
countryFinland   0.491725   0.081517   6.032155   0e+00
countryNorway    0.385830   0.094254   4.093497   0e+00
countrySweden    0.433789   0.060648   7.152599   0e+00
atanh(rho) MZ    0.884166   0.099366   8.898113   0e+00
atanh(rho) DZ    0.271770   0.073240   3.710668   2e-04

 MZ/DZ     Complete pairs MZ/DZ
 4716/8835 997/1809            

                           Estimate 2.5%    97.5%  
Tetrachoric correlation MZ 0.70850  0.59760 0.79280
Tetrachoric correlation DZ 0.26527  0.12752 0.39298

MZ:
                     Estimate 2.5%    97.5%  
Concordance          0.02347  0.01664 0.03300
Casewise Concordance 0.41255  0.32395 0.50721
Marginal             0.05688  0.04643 0.06953
Rel.Recur.Risk       7.25251  5.40099 9.10403
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.00794  0.00489 0.01287
Casewise Concordance 0.13966  0.09312 0.20421
Marginal             0.05688  0.04643 0.06953
Rel.Recur.Risk       2.45511  1.47912 3.43110

                         Estimate 2.5%    97.5%  
Broad-sense heritability 0.88646  0.55608 1.21683

Stratified analysis

bp8 <- bptwin(cancer~strata(country),zyg="zyg",DZ="DZ",id="id",
              type="u",data=prtw,weight="w")
summary(bp8)
------------------------------------------------------------
Strata 'Denmark'

                Estimate    Std.Err          Z p-value
(Intercept)    -1.583608   0.051241 -30.904856  0.0000
atanh(rho) MZ   0.992896   0.217349   4.568215  0.0000
atanh(rho) DZ   0.070588   0.186956   0.377566  0.7058

 MZ/DZ     Complete pairs MZ/DZ
 1334/2789 287/589             

                           Estimate 2.5%     97.5%   
Tetrachoric correlation MZ  0.75859  0.51308  0.88937
Tetrachoric correlation DZ  0.07047 -0.28750  0.41117

MZ:
                     Estimate 2.5%     97.5%   
Concordance           0.02611  0.01584  0.04274
Casewise Concordance  0.46093  0.28426  0.64799
Marginal              0.05664  0.04623  0.06922
Rel.Recur.Risk        8.13766  4.72047 11.55486
DZ:
                     Estimate 2.5%     97.5%   
Concordance           0.00420  0.00110  0.01596
Casewise Concordance  0.07422  0.01888  0.25037
Marginal              0.05664  0.04623  0.06922
Rel.Recur.Risk        1.31043 -0.43515  3.05601

                         Estimate 2.5% 97.5%
Broad-sense heritability   1      NaN  NaN  

------------------------------------------------------------
Strata 'Finland'

                Estimate    Std.Err          Z p-value
(Intercept)    -1.087902   0.063221 -17.207912  0.0000
atanh(rho) MZ   0.859335   0.302752   2.838410  0.0045
atanh(rho) DZ   0.393145   0.179942   2.184840  0.0289

 MZ/DZ    Complete pairs MZ/DZ
 660/1633 134/316             

                           Estimate 2.5%    97.5%  
Tetrachoric correlation MZ 0.69592  0.25985 0.89623
Tetrachoric correlation DZ 0.37407  0.04044 0.63265

MZ:
                     Estimate 2.5%    97.5%  
Concordance          0.07008  0.03975 0.12064
Casewise Concordance 0.50666  0.27641 0.73412
Marginal             0.13832  0.11316 0.16801
Rel.Recur.Risk       3.66298  1.85349 5.47246
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.04160  0.02237 0.07607
Casewise Concordance 0.30073  0.16558 0.48242
Marginal             0.13832  0.11316 0.16801
Rel.Recur.Risk       2.17417  1.00995 3.33838

                         Estimate 2.5%     97.5%   
Broad-sense heritability  0.64369 -0.21675  1.50414

------------------------------------------------------------
Strata 'Norway'

                Estimate    Std.Err          Z p-value
(Intercept)    -1.192293   0.079124 -15.068598  0.0000
atanh(rho) MZ   0.916471   0.301133   3.043409  0.0023
atanh(rho) DZ   0.533761   0.252070   2.117509  0.0342

 MZ/DZ   Complete pairs MZ/DZ
 617/928 115/155             

                           Estimate 2.5%    97.5%  
Tetrachoric correlation MZ 0.72422  0.31516 0.90635
Tetrachoric correlation DZ 0.48825  0.03969 0.77303

MZ:
                     Estimate 2.5%    97.5%  
Concordance          0.05918  0.03218 0.10633
Casewise Concordance 0.50764  0.27633 0.73572
Marginal             0.11657  0.08945 0.15057
Rel.Recur.Risk       4.35466  2.15709 6.55223
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.03945  0.01840 0.08257
Casewise Concordance 0.33842  0.15583 0.58636
Marginal             0.11657  0.08945 0.15057
Rel.Recur.Risk       2.90310  0.89710 4.90911

                         Estimate 2.5%     97.5%   
Broad-sense heritability  0.47195 -0.47133  1.41522

------------------------------------------------------------
Strata 'Sweden'

                Estimate    Std.Err          Z p-value
(Intercept)    -1.149412   0.032155 -35.745836  0.0000
atanh(rho) MZ   0.836864   0.125476   6.669520  0.0000
atanh(rho) DZ   0.199677   0.092907   2.149202  0.0316

 MZ/DZ     Complete pairs MZ/DZ
 2105/3485 461/749             

                           Estimate 2.5%    97.5%  
Tetrachoric correlation MZ 0.68414  0.53057 0.79423
Tetrachoric correlation DZ 0.19706  0.01758 0.36425

MZ:
                     Estimate 2.5%    97.5%  
Concordance          0.06055  0.04659 0.07835
Casewise Concordance 0.48365  0.38001 0.58872
Marginal             0.12519  0.11277 0.13877
Rel.Recur.Risk       3.86327  3.00137 4.72517
DZ:
                     Estimate 2.5%    97.5%  
Concordance          0.02515  0.01672 0.03766
Casewise Concordance 0.20088  0.13541 0.28746
Marginal             0.12519  0.11277 0.13877
Rel.Recur.Risk       1.60452  0.99901 2.21004

                         Estimate 2.5%    97.5%  
Broad-sense heritability 0.97416  0.53594 1.41238

Wald test (stratified vs main effect)

B <- contr(3,4)[-(1:3),]
compare(bp8,contrast=B)
	- Wald test -

	Null Hypothesis:
	[Denmark.atanh(rho) MZ] - [Finland.atanh(rho) MZ] = 0
	[Denmark.atanh(rho) MZ] - [Norway.atanh(rho) MZ] = 0
	[Denmark.atanh(rho) MZ] - [Sweden.atanh(rho) MZ] = 0
	[Denmark.atanh(rho) DZ] - [Finland.atanh(rho) DZ] = 0
	[Denmark.atanh(rho) DZ] - [Norway.atanh(rho) DZ] = 0
	[Denmark.atanh(rho) DZ] - [Sweden.atanh(rho) DZ] = 0

data:  
chisq = 3.4972, df = 6, p-value = 0.7443
sample estimates:
    Estimate   Std.Err       2.5%     97.5%
  0.13356053 0.3726923 -0.5969029 0.8640239
  0.07642511 0.3713780 -0.6514624 0.8043126
  0.15603178 0.2509676 -0.3358556 0.6479191
 -0.32255628 0.2594839 -0.8311353 0.1860227
 -0.46317298 0.3138347 -1.0782776 0.1519316
 -0.12908846 0.2087690 -0.5382682 0.2800912

Created: 2014-05-09 Fri 12:09