Mixture models

Sometimes we can postulate that the subjects are sampled from a population comprised of multiple subpopulations, each with a distinct set of PK parameter typical values, and there are no covariates data available to indicate which subpopulation a subject is from. For such a scenario, mixture models can be used to estimate $P(m)$, the probability that a subject belongs to subpopulation $m$. In other words, unknown parameter $P(m)$ reflects the proportion of subpopulation $m$ in the entire population. That is why $P$ is known as mixing proportions. In general, with $M$ subpopulations, $$ \sum_{m=1}^M P(m) = 1 $$ Namely, $P$ is on a $(M-1)$-simplex. The individual likelihood $L_i$ is then replaced by $$ L_i=\sum_{m=1}^M p(y_i|\Sigma_i, \theta^{(m)}, \Omega)P(m). $$ Here $\theta^{(m)}$ is $\theta$ for subpopulation $m$. The above $P$-weighted average of the invidual likelihood (over all subpopulations) is known as the mixture likelihood. .

For example, assume the body weight in the previous PopPK problem is not measured, we postulate that the population comprised of two subpopulations (possibly due to some unknown body weight threshold), denoted by superscript $(1)$ and $(2)$: $$ \theta^{(1)} = (\theta_1, \theta_3),\quad \theta^{(2)} = (\theta_1\theta_2, \theta_3\theta_4), $$ so that $$ \text{CL}^{(1)}=\theta_1e^{\eta_1}, \quad \text{Vd}^{(1)}=\theta_3e^{\eta_2}, \quad \text{CL}^{(2)}=\theta_2\theta_1e^{\eta_3}, \quad \text{Vd}^{(2)}=\theta_4\theta_3e^{\eta_4}. $$

We can use $MIX record to describe the model:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30


  $PROB PHENOBARB POPULATION DATA MIXTURE MODEL
  $DATA DATA2 (6F7.0)
  $INPUT ID TIME AMT WT APGR CP=DV
  $SUBROUTINES ADVAN1 TRANS2
  $MIX
  NSPOP=2
  P(1)=THETA(5)
  P(2)=1.-THETA(5)
  $PK
  CALLFL=1
  EST=MIXEST
  CL1=THETA(1)*EXP(ETA(1))
  V1=THETA(3)*EXP(ETA(2))
  CL2=THETA(2)*THETA(1)*EXP(ETA(3))
  V2=THETA(4)*THETA(3)*EXP(ETA(4))
  IF (MIXNUM.EQ.1) THEN
    CL=CL1
    V=V1
  ELSE
    CL=CL2
    V=V2
  ENDIF
  S1=V
  $ERROR
  Y=F*EXP(EPS(1))
  $THETA (0,.0047) (0,2.0) (0,.99) (0,2.0) (0,.5,1)
  $OMEGA BLOCK(2) .05 .01 .03
  $OMEGA BLOCK(2) SAME
  $SIGMA .02
  $ESTM MAXEVAL=500 PRINT=5 INTERACTION METHOD=1

NSPOP=2 in the $MIX block indicates there are two subpopulations ($K=2$). Reserved variable P(m), m=1,2 is for $P(m)$, modeled by a THETA element. User can use the reserved variable MIXNUM as subpopulation index to selectively specify CL and V. By fitting the above model we estimate both sets of typical values and the mixing fraction $P$.

A slightly more succinct way to write mixture models is to use indicator variables. Take the following control stream for example (examples/example5, see also examples/example3 and examples/example4).

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51


;# Model Desc: Population Mixture Problem in 1 Compartment model, with rate constant parameter mean modeled for two sub-populations, but its inter-subject variance is the same in both sub-populations. See
;#  https://nmusers.github.io/docs/reference-manual/control-records/mix/
;# Project Name: nm7examples

$PROB RUN example5 (from ad1tr1m4t)
$INPUT C SET ID JID TIME CONC=DV DOSE=AMT RATE EVID MDV CMT VC1 K101 VC2 K102 SIGZ PROB
$DATA example5.csv IGNORE=C

$SUBROUTINES ADVAN1 TRANS1

$MIX
P(1)=THETA(4)
P(2)=1.0-THETA(4)
NSPOP=2

$PK
Q=1
IF(MIXNUM.EQ.2) Q=0
MU_1=THETA(1)
;# Note that MU_2 can be modeled as THETA(2) or THETA(3), depending on the MIXNUM value.
;# Also, we are avoiding IF/THEN blocks.
MU_2=Q*THETA(2)+(1.0-Q)*THETA(3)
V=DEXP(MU_1+ETA(1))
K=DEXP(MU_2+ETA(2))
S1=V

$ERROR
Y = F + F*EPS(1)

$THETA
(-1000.0  4.3 1000.0)  ;[MU_1]
(-1000.0 -2.9 1000.0)  ;[MU_2-1]
(-1000.0 -0.67 1000.0) ;[MU_2-2]
(0.0001 0.667 0.9999)  ;[P(1)]

$OMEGA BLOCK(2)
0.04 ;[p]
0.01  ;[f]
0.04 ;[p]

$SIGMA
0.01 ;[p]

$EST METHOD=ITS INTERACTION NITER=100 PRINT=1 NOABORT SIGL=8 FILE=example5.ext CTYPE=3
$EST METHOD=IMPMAP INTERACTION NITER=20 ISAMPLE=300 PRINT=1 NOABORT SIGL=8
$EST METHOD=IMP INTERACTION NITER=20 ISAMPLE=1000 PRINT=1 NOABORT SIGL=6
$EST NBURN=500 NITER=500 METHOD=SAEM INTERACTION PRINT=10 SIGL=6 ISAMPLE=2
$EST METHOD=IMP INTERACTION NITER=5 ISAMPLE=1000 PRINT=1 NOABORT SIGL=6 EONLY=1 MAPITER=0
$EST METHOD=BAYES INTERACTION NBURN=2000 NITER=5000 PRINT=10  FILE=example5.txt SIGL=8
$EST MAXEVAL=9999 NSIG=2 SIGL=8 PRINT=10 FILE=example5.ext METHOD=CONDITIONAL INTERACTION NOABORT
$COV MATRIX=R

In this model, indicator variable Q is set according to MIXNUM, and used to define the typical value of $K$.