MixedModels.jl confidence intervals and repeated measure, how to get?

Specific Domains Statistics
1

I have a data:

julia> df
298×6 DataFrame
│ Row │ Subject      │ Period       │ Seq          │ Formulation  │ logData │ Data    │
│     │ Categorical… │ Categorical… │ Categorical… │ Categorical… │ Float64 │ Float64 │
├─────┼──────────────┼──────────────┼──────────────┼──────────────┼─────────┼─────────┤
│ 1   │ 1            │ 1            │ BABA         │ R            │ 7.73454 │ 2285.96 │
│ 2   │ 1            │ 2            │ BABA         │ T            │ 7.57857 │ 1955.82 │
│ 3   │ 1            │ 3            │ BABA         │ R            │ 7.20485 │ 1345.94 │
│ 4   │ 1            │ 4            │ BABA         │ T            │ 7.95726 │ 2856.24 │
⋮
│ 294 │ 77           │ 4            │ BABA         │ T            │ 7.76337 │ 2352.82 │
│ 295 │ 78           │ 1            │ BABA         │ R            │ 9.92199 │ 20373.5 │
│ 296 │ 78           │ 2            │ BABA         │ T            │ 9.52439 │ 13689.6 │
│ 297 │ 78           │ 3            │ BABA         │ R            │ 9.93232 │ 20585.0 │
│ 298 │ 78           │ 4            │ BABA         │ T            │ 10.1064 │ 24498.1 │

I try to use MixedModels.jl:

lm = LinearMixedModel(@formula(logData ~ 1 + Period + Seq + Formulation  + (0 +  Formulation|Subject)), df)
fm1 = fit!(lm, REML=true)

How can i get confidence interval for formulation difference (contrast)?

How can i fit repeated measure and use Period as repeated factor?

2

Not sure I understand correctly but it looks as though there’s only two level to your variable so the reported coefficient estimate (and standard error) is the difference between T and R?

3

Hello!

Yes! I try to make CI for T-R difference.

4

Exactly I want to get intra-subject variance for T and R separately. But for code above i get results:

inear mixed model fit by REML
 Formula: logData ~ 1 + Period + Seq + Formulation + ((0 + Formulation) | Subject)
 REML criterion at convergence: 535.1437908228327

Variance components:
               Column      Variance   Std.Dev.    Corr.
 Subject  Formulation: R  0.75000775 0.86602988
          Formulation: T  0.66580315 0.81596762  1.00
 Residual                 0.15924154 0.39905079
 Number of obs: 298; levels of grouping factors: 77

Much higher than in SAS or SPSS. I suppose it is wrong fitting because repeated effect not included.

May be this fitting more correct:

lm = LinearMixedModel(@formula(logData ~ 1 + Period + Seq + Formulation  +  (0 + Formulation|Subject) + (0 + Formulation|Period)), df)

But this results far from SPSS:
berep

5

I note that the difference between the coefficients is the same at 0.085 for Julia and SPSS - might be due to different selection of base category for the categorical coding? Have a look at Contrast coding categorical variables · StatsModels.jl

Are you able to share your data?

6

Yep!

Subject,Period,Seq,Formulation,logData,Data
1,1,BABA,R,7.734541,2285.96
1,2,BABA,T,7.578565,1955.82
1,3,BABA,R,7.204848,1345.94
1,4,BABA,T,7.957261,2856.24
2,1,ABAB,T,8.055704,3151.72
2,2,ABAB,R,7.859143,2589.3
2,3,ABAB,T,8.004011,2992.94
2,4,ABAB,R,7.788792,2413.4
3,1,ABAB,T,8.090935,3264.74
3,2,ABAB,R,8.088844,3257.92
3,3,ABAB,T,8.039332,3100.54
3,4,ABAB,R,8.037272,3094.16
4,1,ABAB,T,7.095363,1206.36
4,2,ABAB,R,7.175153,1306.56
4,3,ABAB,T,7.367153,1583.12
4,4,ABAB,R,7.207445,1349.44
5,1,BABA,R,8.263822,3880.9
5,2,BABA,T,8.898759,7322.88
5,3,BABA,R,8.396078,4429.66
5,4,BABA,T,8.108587,3322.88
6,1,BABA,R,6.885591,978.08
6,2,BABA,T,7.099235,1211.04
6,3,BABA,R,6.881288,973.88
6,4,BABA,T,7.048213,1150.8
7,1,ABAB,T,7.980728,2924.06
7,2,ABAB,R,7.736298,2289.98
7,3,ABAB,T,7.821755,2494.28
7,4,ABAB,R,8.083063,3239.14
8,1,BABA,R,7.793776,2425.46
8,2,BABA,T,8.217638,3705.74
8,3,BABA,R,7.544893,1891.06
8,4,BABA,T,9.102657,8979.12
9,1,BABA,R,8.249319,3825.02
9,2,BABA,T,8.578296,5315.04
9,3,BABA,R,8.66788,5813.16
9,4,BABA,T,9.348004,11475.9
10,1,BABA,R,8.321728,4112.26
10,2,BABA,T,8.248754,3822.86
10,3,BABA,R,7.807843,2459.82
10,4,BABA,T,8.437448,4616.76
11,1,ABAB,T,8.061581,3170.3
11,2,ABAB,R,7.726257,2267.1
11,4,ABAB,R,7.440335,1703.32
12,1,ABAB,T,8.005427,2997.18
12,2,ABAB,R,7.991179,2954.78
12,3,ABAB,T,8.56649,5252.66
12,4,ABAB,R,8.228054,3744.54
13,1,ABAB,T,7.628372,2055.7
13,2,ABAB,R,6.890914,983.3
13,3,ABAB,T,7.479469,1771.3
13,4,ABAB,R,8.099609,3293.18
14,1,BABA,R,7.371879,1590.62
14,2,BABA,T,7.040134,1141.54
14,3,BABA,R,7.121527,1238.34
14,4,BABA,T,7.159136,1285.8
15,1,ABAB,T,7.293358,1470.5
15,2,ABAB,R,7.023616,1122.84
15,3,ABAB,T,7.372859,1592.18
15,4,ABAB,R,7.469175,1753.16
16,1,BABA,R,7.542288,1886.14
16,2,BABA,T,7.638815,2077.28
16,3,BABA,R,7.69513,2197.62
16,4,BABA,T,7.693773,2194.64
17,1,ABAB,T,6.444386,629.16
17,2,ABAB,R,6.211283,498.34
17,3,ABAB,T,6.313077,551.74
17,4,ABAB,R,5.945892,382.18
18,1,BABA,R,6.141951,464.96
18,2,BABA,T,7.989506,2949.84
18,3,BABA,R,7.094716,1205.58
18,4,BABA,T,7.671342,2145.96
19,1,BABA,R,7.54394,1889.26
19,2,BABA,T,8.671996,5837.14
19,3,BABA,R,7.773106,2375.84
19,4,BABA,T,7.422649,1673.46
20,1,ABAB,T,6.676378,793.44
20,2,ABAB,R,7.06452,1169.72
20,4,ABAB,R,6.978027,1072.8
21,1,BABA,R,7.642898,2085.78
21,2,BABA,T,7.771995,2373.2
21,3,BABA,R,7.350516,1557
21,4,BABA,T,7.666353,2135.28
22,1,BABA,R,5.663169,288.06
22,2,BABA,T,5.736508,309.98
22,3,BABA,R,5.781299,324.18
22,4,BABA,T,5.728735,307.58
23,1,ABAB,T,6.263017,524.8
23,2,ABAB,R,5.921149,372.84
23,3,ABAB,T,6.25175,518.92
23,4,ABAB,R,6.404501,604.56
24,1,ABAB,T,8.677088,5866.94
24,3,ABAB,T,8.621153,5547.78
24,4,ABAB,R,8.386355,4386.8
25,1,ABAB,T,8.296162,4008.46
25,2,ABAB,R,7.548998,1898.84
25,3,ABAB,T,7.355782,1565.22
25,4,ABAB,R,8.491941,4875.32
26,1,ABAB,T,7.087958,1197.46
26,2,ABAB,R,5.801574,330.82
26,3,ABAB,T,7.151611,1276.16
26,4,ABAB,R,5.97843,394.82
27,1,BABA,R,9.534102,13823.18
27,2,BABA,T,8.938377,7618.82
27,3,BABA,R,9.158346,9493.34
27,4,BABA,T,9.096997,8928.44
28,1,BABA,R,6.846794,940.86
28,2,BABA,T,7.080616,1188.7
28,3,BABA,R,6.782215,882.02
28,4,BABA,T,7.111822,1226.38
29,1,BABA,R,7.684894,2175.24
29,2,BABA,T,7.883959,2654.36
29,3,BABA,R,8.081865,3235.26
29,4,BABA,T,8.017406,3033.3
30,1,ABAB,T,7.085818,1194.9
30,2,ABAB,R,6.717393,826.66
30,3,ABAB,T,6.414082,610.38
30,4,ABAB,R,6.387115,594.14
31,1,BABA,R,8.320857,4108.68
31,2,BABA,T,8.90917,7399.52
31,4,BABA,T,8.403267,4461.62
32,1,ABAB,T,6.674839,792.22
32,2,ABAB,R,6.907495,999.74
32,3,ABAB,T,7.072761,1179.4
32,4,ABAB,R,7.42593,1678.96
33,1,BABA,R,8.275254,3925.52
33,2,BABA,T,8.240053,3789.74
33,3,BABA,R,8.150127,3463.82
33,4,BABA,T,8.42872,4576.64
34,1,BABA,R,7.443418,1708.58
34,2,BABA,T,7.824382,2500.84
34,3,BABA,R,7.141483,1263.3
34,4,BABA,T,7.624824,2048.42
35,1,ABAB,T,6.84913,943.06
35,2,ABAB,R,6.645377,769.22
35,3,ABAB,T,6.743824,848.8
35,4,ABAB,R,7.084964,1193.88
36,1,ABAB,T,7.840085,2540.42
36,2,ABAB,R,7.645484,2091.18
36,3,ABAB,T,7.856962,2583.66
36,4,ABAB,R,7.597888,1993.98
37,1,ABAB,T,6.746929,851.44
37,2,ABAB,R,6.482924,653.88
37,3,ABAB,T,7.771194,2371.3
37,4,ABAB,R,7.150999,1275.38
38,1,BABA,R,8.7086,6054.76
38,2,BABA,T,8.898663,7322.18
38,3,BABA,R,8.81685,6746.98
38,4,BABA,T,8.872165,7130.7
39,1,ABAB,T,8.670024,5825.64
39,2,ABAB,R,8.773821,6462.82
39,3,ABAB,T,8.9093,7400.48
39,4,ABAB,R,8.731795,6196.84
40,1,BABA,R,7.432732,1690.42
40,2,BABA,T,7.164643,1292.9
40,3,BABA,R,7.328043,1522.4
40,4,BABA,T,6.972213,1066.58
41,1,BABA,R,7.931306,2783.06
41,2,BABA,T,7.046717,1149.08
41,3,BABA,R,6.777555,877.92
41,4,BABA,T,6.349873,572.42
42,1,ABAB,T,8.467805,4759.06
42,2,ABAB,R,8.671102,5831.92
42,4,ABAB,R,8.33201,4154.76
43,1,ABAB,T,8.594021,5399.28
43,2,ABAB,R,8.598939,5425.9
43,3,ABAB,T,8.376666,4344.5
43,4,ABAB,R,8.413396,4507.04
44,1,ABAB,T,8.632502,5611.1
44,2,ABAB,R,8.602295,5444.14
44,3,ABAB,T,8.4776,4805.9
44,4,ABAB,R,8.509294,4960.66
45,1,BABA,R,6.561992,707.68
45,2,BABA,T,8.211119,3681.66
45,3,BABA,R,9.823051,18454.26
45,4,BABA,T,6.911209,1003.46
46,1,ABAB,T,7.783491,2400.64
46,2,ABAB,R,7.258835,1420.6
46,3,ABAB,T,7.044626,1146.68
46,4,ABAB,R,8.518337,5005.72
47,1,BABA,R,6.180182,483.08
47,2,BABA,T,6.940513,1033.3
47,3,BABA,R,6.468537,644.54
47,4,BABA,T,6.515157,675.3
48,1,BABA,R,7.676511,2157.08
48,2,BABA,T,8.044742,3117.36
48,3,BABA,R,7.943122,2816.14
48,4,BABA,T,7.955215,2850.4
49,1,ABAB,T,9.565322,14261.54
49,2,ABAB,R,10.184506,26489.56
49,3,ABAB,T,10.065847,23525.66
49,4,ABAB,R,9.963818,21243.76
50,1,ABAB,T,7.347454,1552.24
50,2,ABAB,R,7.358398,1569.32
50,3,ABAB,T,7.644919,2090
50,4,ABAB,R,7.299784,1479.98
51,1,BABA,R,8.251779,3834.44
51,2,BABA,T,8.496942,4899.76
51,3,BABA,R,8.216872,3702.9
51,4,BABA,T,8.644182,5677.02
52,1,BABA,R,8.687091,5925.92
52,2,BABA,T,6.875129,967.9
52,3,BABA,R,6.68088,797.02
52,4,BABA,T,6.84522,939.38
53,1,ABAB,T,8.168622,3528.48
53,2,ABAB,R,7.61941,2037.36
53,3,ABAB,T,8.074549,3211.68
53,4,ABAB,R,7.974787,2906.74
54,1,BABA,R,6.842854,937.16
54,2,BABA,T,8.752733,6327.96
54,3,BABA,R,6.96122,1054.92
54,4,BABA,T,7.476484,1766.02
55,1,ABAB,T,8.142639,3437.98
55,2,ABAB,R,8.224646,3731.8
55,3,ABAB,T,8.483165,4832.72
55,4,ABAB,R,8.104776,3310.24
56,1,ABAB,T,6.918834,1011.14
56,2,ABAB,R,6.483138,654.02
56,3,ABAB,T,6.75528,858.58
56,4,ABAB,R,6.811377,908.12
57,1,BABA,R,6.91109,1003.34
57,2,BABA,T,8.46378,4739.94
57,3,BABA,R,6.54799,697.84
57,4,BABA,T,7.825852,2504.52
58,1,BABA,R,8.778994,6496.34
58,2,BABA,T,8.691039,5949.36
58,3,BABA,R,8.700078,6003.38
58,4,BABA,T,8.759939,6373.72
59,1,BABA,R,7.128961,1247.58
59,2,BABA,T,7.018294,1116.88
59,3,BABA,R,7.061969,1166.74
59,4,BABA,T,7.885472,2658.38
60,1,ABAB,T,10.432044,33929.62
60,2,ABAB,R,10.124366,24943.44
60,3,ABAB,T,9.857979,19110.22
60,4,ABAB,R,9.457605,12805.18
62,1,ABAB,T,7.73215,2280.5
62,2,ABAB,R,7.446865,1714.48
62,3,ABAB,T,8.302583,4034.28
62,4,ABAB,R,8.137618,3420.76
63,1,ABAB,T,8.12466,3376.72
63,2,ABAB,R,7.71548,2242.8
63,3,ABAB,T,7.449812,1719.54
63,4,ABAB,R,7.758897,2342.32
64,1,BABA,R,6.816013,912.34
64,2,BABA,T,7.651795,2104.42
64,3,BABA,R,7.630966,2061.04
64,4,BABA,T,7.310884,1496.5
65,1,BABA,R,8.283479,3957.94
65,2,BABA,T,8.681962,5895.6
65,3,BABA,R,8.675833,5859.58
65,4,BABA,T,8.531782,5073.48
66,1,ABAB,T,7.061077,1165.7
66,2,ABAB,R,7.129794,1248.62
66,3,ABAB,T,7.06363,1168.68
66,4,ABAB,R,7.170443,1300.42
67,1,BABA,R,7.087908,1197.4
67,2,BABA,T,7.020495,1119.34
68,1,BABA,R,7.444085,1709.72
68,2,BABA,T,7.836923,2532.4
68,3,BABA,R,7.365825,1581.02
68,4,BABA,T,7.940014,2807.4
69,1,ABAB,T,7.93696,2798.84
69,2,ABAB,R,7.805515,2454.1
69,4,ABAB,R,8.582014,5334.84
70,1,BABA,R,8.370645,4318.42
70,2,BABA,T,7.6883,2182.66
70,3,BABA,R,7.408021,1649.16
70,4,BABA,T,7.390379,1620.32
71,1,ABAB,T,6.153243,470.24
71,2,ABAB,R,5.33773,208.04
72,1,ABAB,T,7.648883,2098.3
72,2,ABAB,R,7.559955,1919.76
72,3,ABAB,T,7.943698,2817.76
72,4,ABAB,R,7.621195,2041
73,1,ABAB,T,8.804973,6667.32
73,2,ABAB,R,8.573543,5289.84
73,3,ABAB,T,8.895668,7300.28
73,4,ABAB,R,9.181101,9711.84
74,1,BABA,R,7.619116,2036.76
74,2,BABA,T,7.574579,1948.04
74,3,BABA,R,7.339265,1539.58
74,4,BABA,T,7.63971,2079.14
75,1,ABAB,T,6.642878,767.3
75,2,ABAB,R,6.953015,1046.3
75,3,ABAB,T,7.237318,1390.36
75,4,ABAB,R,8.012741,3019.18
76,1,ABAB,T,9.400754,12097.5
76,2,ABAB,R,9.448918,12694.42
76,3,ABAB,T,9.305581,10999.24
76,4,ABAB,R,9.149158,9406.52
77,1,BABA,R,7.017058,1115.5
77,2,BABA,T,7.016879,1115.3
77,3,BABA,R,7.013718,1111.78
77,4,BABA,T,7.76337,2352.82
78,1,BABA,R,9.921992,20373.54
78,2,BABA,T,9.524392,13689.6
78,3,BABA,R,9.932319,20585.02
78,4,BABA,T,10.106352,24498.14
1 Like
7

Notice the estimate of the within-Subject correlation of the random effects from MixedModels is 1.00. This is an indication that the estimate of the covariance matrix is singular.

I don’t know how SAS and SPSS formulate the model and estimate the parameters because their code is proprietary. Often, and I suspect it is the case here, the inverse of the covariance matrix is used in the calculation. Naturally this becomes problematic when trying to evaluate at a singular covariance matrix.

The methods in the lme4 package for R and the MixedModels package for Julia were specifically designed to avoid this problem

1 Like
8

You probably don’t want to try to fit a model with a term like (0 + Formulation|Period) in it because there are only 4 levels of Period. You can’t expect to estimate two variances and a covariance from 4 distinct levels of Period.

2 Likes
9

Douglas, thank you so much for reply! I understand that term (0 + Formulation|Period) not appropriate in this situation :frowning: I try to fit (0 + Formulation|Subject*Period) - it leads to error. So I make Subject*Period manualy and call is SP:

Linear mixed model fit by REML
 Formula: logData ~ 1 + Period + Seq + Formulation + ((0 + Formulation) | Subject) + ((0 + Formulation) | SP)
 REML criterion at convergence: 530.1445130337493

Variance components:
               Column       Variance    Std.Dev.    Corr.
 SP       Formulation: R  0.0903490313 0.30058116
          Formulation: T  0.0056246077 0.07499738  0.95
 Subject  Formulation: R  0.7276018429 0.85299580
          Formulation: T  0.6862586172 0.82840728  1.00
 Residual                 0.1117692954 0.33431915
 Number of obs: 298; levels of grouping factors: 298, 77

  Fixed-effects parameters:
─────────────────────────────────────────────────────────
                  Estimate  Std.Error    z value  P(>|z|)
─────────────────────────────────────────────────────────
(Intercept)     7.62921     0.145582   52.4047     <1e-99
Period: 2       0.0446378   0.0636371   0.701443   0.4830
Period: 3       0.00659022  0.0644218   0.102298   0.9185
Period: 4       0.0727199   0.0639714   1.13676    0.2556
Seq: BABA       0.0204706   0.197223    0.103794   0.9173
Formulation: T  0.145464    0.0465012   3.12818    0.0018
─────────────────────────────────────────────────────────

maybe this fitting wrong to…

For SPSS code looke like this:

MIXED logData BY Period Formulation Seq Subject
  /CRITERIA=CIN(90) MXITER(200) MXSTEP(20) SCORING(2) SINGULAR(0.000000000001) HCONVERGE(0,
    RELATIVE) LCONVERGE(0.0000000000001, RELATIVE) PCONVERGE(0, RELATIVE)
  /FIXED=Period Formulation Sequence | SSTYPE(3)
  /METHOD=REML
  /RANDOM=Formulation | SUBJECT(Subject) COVTYPE(CSH)
  /REPEATED=Formulation | SUBJECT(Subject*Period) COVTYPE(DIAG)
  /EMMEANS=TABLES(Formulation) COMPARE REFCAT(FIRST) ADJ(LSD).

Ok. So. Is it possible to get CI for difference T-R from MixedModels?

10

I’m sorry but I haven’t used either SPSS or SAS in about 40 years and never did understand the REPEATED formulation in the first place. You will need to get someone who understands what the SPSS model means.

1 Like