You are here

Parameter contraints across subgroups

4 posts / 0 new
Last post
Julia's picture
Offline
Joined: 03/29/2012 - 09:40
Parameter contraints across subgroups

Hi all,

I am running a simple univariate saturated threshold model to find out the effects of covariates on the variable of interest across three groups - MZ twins, DZ twins and unrelated individuals. The regressions are the following: y ~age + sex + IBS.
As one of the research questions, I would like to specify such a constraint for IBS coefficients: beta_unrel = 2*beta_DZ. This is where I got stuck. I looked at the following thread https://openmx.ssri.psu.edu/thread/446 and wanted to specify an mxAlgebra instead of mxConstraint, but I struggle to modify my original model to test this assumption.

This is the chunk of the code for the regressions from the original model:

meanMZ <- mxMatrix( type="Zero", nrow=1, ncol=ntv, name="expMeanMZ" )
meanDZ <- mxMatrix( type="Zero", nrow=1, ncol=ntv, name="expMeanDZ" )
meanUnrel <- mxMatrix( type="Zero", nrow=1, ncol=ntv, name="expMeanUnrel" )
 
threMZ  <- mxMatrix( type="Full", nrow=1, ncol=ntv, free=TRUE, values=.8, labels='threMZ', name="threMZ" )
threDZ  <- mxMatrix( type="Full", nrow=1, ncol=ntv, free=TRUE, values=.8, labels='threDZ', name="threDZ" )
threUnrel   <- mxMatrix( type="Full", nrow=1, ncol=ntv, free=TRUE, values=.8, labels='threUnrel', name="threUnrel" )
 
defAge1     <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels="data.Zage_s_1", name="Age1" )
defAge2     <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels="data.Zage_s_2", name="Age2" )
 
B_AgeMZ         <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaAgeMZ", name="bAgeMZ" )
B_AgeDZ         <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaAgeDZ", name="bAgeDZ" )
B_AgeUnrel      <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaAgeUnrel", name="bAgeUnrel" )
 
threAge1MZ   <- mxAlgebra(  Age1%*%bAgeMZ, name="AgeR1MZ")
threAge2MZ   <- mxAlgebra(  Age2%*%bAgeMZ, name="AgeR2MZ")
threAge1DZ   <- mxAlgebra(  Age1%*%bAgeDZ, name="AgeR1DZ")
threAge2DZ   <- mxAlgebra(  Age2%*%bAgeDZ, name="AgeR2DZ")
threAge1Unrel   <- mxAlgebra(  Age1%*%bAgeUnrel, name="AgeR1Unrel")
threAge2Unrel   <- mxAlgebra(  Age2%*%bAgeUnrel, name="AgeR2Unrel")
 
defSex1     <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels="data.sex_s_1", name="Sex1" )
defSex2     <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels="data.sex_s_2", name="Sex2" )
 
B_SexMZ         <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaSexMZ", name="bSexMZ" )
B_SexDZ         <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaSexDZ", name="bSexDZ" )
B_SexUnrel      <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaSexUnrel", name="bSexUnrel" )
 
threSex1MZ   <- mxAlgebra(  Sex1%*%bSexMZ, name="SexR1MZ")
threSex2MZ   <- mxAlgebra(  Sex2%*%bSexMZ, name="SexR2MZ")
threSex1DZ   <- mxAlgebra(  Sex1%*%bSexDZ, name="SexR1DZ")
threSex2DZ   <- mxAlgebra(  Sex2%*%bSexDZ, name="SexR2DZ")
threSex1Unrel   <- mxAlgebra(  Sex1%*%bSexUnrel, name="SexR1Unrel")
threSex2Unrel   <- mxAlgebra(  Sex2%*%bSexUnrel, name="SexR2Unrel")
 
 
defIBS1     <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels="data.ibs_gen_s_1", name="IBS1" )
defIBS2     <- mxMatrix( type="Full", nrow=1, ncol=1, free=FALSE, labels="data.ibs_gen_s_2", name="IBS2" )
 
B_IBSMZ         <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaIBSMZ", name="bIBSMZ" )
B_IBSDZ         <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaIBSDZ", name="bIBSDZ" )
B_IBSUnrel      <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaIBSUnrel", name="bIBSUnrel" )
 
threIBS1MZ   <- mxAlgebra(  IBS1%*%bIBSMZ, name="IBSR1MZ")
threIBS2MZ   <- mxAlgebra(  IBS2%*%bIBSMZ, name="IBSR2MZ")
threIBS1DZ   <- mxAlgebra(  IBS1%*%bIBSDZ, name="IBSR1DZ")
threIBS2DZ   <- mxAlgebra(  IBS2%*%bIBSDZ, name="IBSR2DZ")
threIBS1Unrel   <- mxAlgebra(  IBS1%*%bIBSUnrel, name="IBSR1Unrel")
threIBS2Unrel   <- mxAlgebra(  IBS2%*%bIBSUnrel, name="IBSR2Unrel")
 
expThreMZ   <- mxAlgebra( expression= threMZ + cbind(AgeR1MZ,AgeR2MZ) + cbind(SexR1MZ,SexR2MZ) + cbind(IBSR1MZ,IBSR2MZ), name="expThreMZ" )
expThreDZ   <- mxAlgebra( expression= threDZ + cbind(AgeR1DZ,AgeR2DZ) + cbind(SexR1DZ,SexR2DZ) + cbind(IBSR1DZ,IBSR2DZ), name="expThreDZ" )
expThreUnrel   <- mxAlgebra( expression= threUnrel + cbind(AgeR1Unrel,AgeR2Unrel) + cbind(SexR1Unrel,SexR2Unrel) + cbind(IBSR1Unrel,IBSR2Unrel), name="expThreUnrel" )

All this, together with other entities, go into three subgroups and into the total model:

modelMZ   <- mxModel( name="MZ", meanMZ, covMZ,
                      threMZ, defAge1, defAge2, B_AgeMZ, threAge1MZ, threAge2MZ, defSex1, defSex2, B_SexMZ, threSex1MZ, threSex2MZ, defIBS1, defIBS2, B_IBSMZ, threIBS1MZ, threIBS2MZ,
                      expThreMZ, dataMZ, expMZ, funML )
modelDZ   <- mxModel( name="DZ", meanDZ, covDZ,
                      threDZ, defAge1, defAge2, B_AgeDZ, threAge1DZ, threAge2DZ, defSex1, defSex2, B_SexDZ, threSex1DZ, threSex2DZ, defIBS1, defIBS2, B_IBSDZ, threIBS1DZ, threIBS2DZ,
                      expThreDZ, dataDZ, expDZ, funML )
modelUnrel   <- mxModel( name="Unrel", meanUnrel, covUnrel,
                      threUnrel, defAge1, defAge2, B_AgeUnrel, threAge1Unrel, threAge2Unrel, defSex1, defSex2, B_SexUnrel, threSex1Unrel, threSex2Unrel, defIBS1, defIBS2, B_IBSUnrel, threIBS1Unrel, threIBS2Unrel,
                      expThreUnrel, dataUnrel, expUnrel, funML )
 
multi     <- mxFitFunctionMultigroup( c("MZ","DZ","Unrel") )
 
univSatModel <- mxModel("univSat", modelMZ, modelDZ, modelUnrel, multi)
univSatFit   <- mxRun(univSatModel)

Now I would like to test whether we can impose betaIBSUnrel = 2*betaIBSDZ. I tried the following:

restrModel = mxModel(model1, name='restrM')
restrModel = omxSetParameters(restrModel , labels='betaIBSDZ', values=0.01, free=T, newlabels = 'betaIBSDZ')
B_IBSDZ_Unrel       <- mxMatrix( type="Full", nrow=1, ncol=nv, free=TRUE, values=.01, label="betaIBSDZ", name="bIBSDZ" )
B_IBSUnrel_matrix   <- mxMatrix( type="Full", nrow=1, ncol=nv, values=.01, label="bIBSUnrel[1,1]", name="bIBSUnrel_matrix" )
B_IBSUnrel_alg = mxAlgebra(2*betaIBSDZ, name='bIBSUnrel')
restrModel = mxModel(GenIBSModel, B_IBSDZ_Unrel, B_IBSUnrel_matrix, B_IBSUnrel_alg)
restrFit = mxRun(restrModel)

But this doesn't seem to do the right thing. The new 'bIBSUnrel' from B_IBSUnrel_alg is calculated along the old 'bIBSUnrel' in the Unrel subgroup. This is the summary of the model fit:

> summary(restrFit)
Summary of restrM
 
free parameters:
           name          matrix             row             col    Estimate  Std.Error A
1     betaIBSDZ          bIBSDZ               1               1 -0.13358435 0.13235843  
2     covT1T2MZ     MZ.expCovMZ strain_dich_s_1 strain_dich_s_2  0.42856769 0.11593230  
3     threshAll       MZ.threMZ               1               1  0.64756220 0.10746847  
4     betaAgeMZ       MZ.bAgeMZ               1               1  0.07123071 0.07590719  
5     betaSexMZ       MZ.bSexMZ               1               1 -0.14514120 0.14911503  
6     covT1T2DZ     DZ.expCovDZ strain_dich_s_1 strain_dich_s_2  0.22069379 0.13846191  
7     betaAgeDZ       DZ.bAgeDZ               1               1  0.37860942 0.09224387  
8     betaSexDZ       DZ.bSexDZ               1               1  0.07751815 0.14411199  
9   threshUnrel Unrel.threUnrel               1               1  0.96161404 0.03937058  
10 betaAgeUnrel Unrel.bAgeUnrel               1               1  0.10529660 0.02691893  
11 betaSexUnrel Unrel.bSexUnrel               1               1 -0.13111884 0.05280346  
12 betaIBSUnrel Unrel.bIBSUnrel               1               1 -0.49658313 0.09702530  
 
> mxEval(c(betaIBSDZ, bIBSUnrel), restrFit)
           [,1]
[1,] -0.1335843
[2,] -0.2671687

You can see that there is still a matrix bIBSUnrel in Unrel group which element was estimated freely (-0.497), whereas the new bIBSUnrel matrix that I created for the submodel is equal to two times betaIBSDZ.

Could you please help me to code this across group constraint and to run a restricted model with it?
Thank you in advance!
Julia

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
namespace

Your syntax doesn't work because B_IBSUnrel_matrix is going into the container model, not the submodel modelUnrel, so the submodel's expectation never "sees" it. Plus, your syntax never fixes the parameter labeled "bIBSUnrel", which is still there in that submodel. I think the easiest thing is probably to just define the regression coefficient for IBS among unrelateds as an algebra,

B_IBSUnrel_alg <- mxAlgebra(2*betaIBSDZ, name='bIBSUnrel')

, or perhaps,

B_IBSUnrel_alg <- mxAlgebra(2*DZ.bIBSDZ, name='bIBSUnrel')

, then redefine modelUnrel to use it instead of B_IBSUnrel,

modelUnrel   <- mxModel( name="Unrel", 
     meanUnrel, covUnrel, threUnrel, defAge1, defAge2, B_AgeUnrel, threAge1Unrel, threAge2Unrel, defSex1, defSex2, 
     B_SexUnrel, threSex1Unrel, threSex2Unrel, defIBS1, defIBS2, B_IBSUnrel_alg, threIBS1Unrel, threIBS2Unrel,
     expThreUnrel, dataUnrel, expUnrel, funML )

and then rebuild the multigroup model using the redefined modelUnrel.

Julia's picture
Offline
Joined: 03/29/2012 - 09:40
Great!

Hi Rob. Thanks a bunch!!! This does seem to work! Thank you again!

AdminRobK's picture
Offline
Joined: 01/24/2014 - 12:15
Glad to hear it!

Glad to hear it!