@@ -119,25 +119,25 @@ static double event_rates
119119 // total pi = (I1-0.5*ell1)/I1
120120 assert (S2 >=0 && I1 >=ell1 && ell1 >=0 );
121121 alpha = (N1 > 0 ) ? Beta21 * S2 * I1 /N1 : 0 ;
122- pi = (I1 > 0 ) ? (I1 - 0.5 * ell1 )/I1 : 0 ;
122+ pi = (I1 > 0 ) ? (I1 - 0.5 * ell1 )/I1 : 1 ;
123123 event_rate += (* rate = alpha * pi ); rate ++ ;
124124 * logpi = log (pi ); logpi ++ ;
125125 assert (R_FINITE (event_rate ));
126126 // 3: Trans_21, s = (0,1)
127127 // s = (0,1): pi = 1/(2*I1), m = ell1
128- pi = 1 - pi ;
128+ pi = ( I1 > 0 ) ? 1 - pi : 1 ;
129129 event_rate += (* rate = alpha * pi ); rate ++ ;
130130 * logpi = log (pi )- log (ell1 ); logpi ++ ;
131131 assert (R_FINITE (event_rate ));
132132 // 4: Trans_12, s = (0,0),(0,1)
133133 assert (S1 >=0 && I2 >=ell2 && ell2 >=0 );
134134 alpha = (N2 > 0 ) ? Beta12 * S1 * I2 /N2 : 0 ;
135- pi = (I2 > 0 ) ? (I2 - ell2 * 0.5 )/I2 : 0 ;
135+ pi = (I2 > 0 ) ? (I2 - ell2 * 0.5 )/I2 : 1 ;
136136 event_rate += (* rate = alpha * pi ); rate ++ ;
137137 * logpi = log (pi ); logpi ++ ;
138138 assert (R_FINITE (event_rate ));
139139 // 5: Trans_12, s = (1,0)
140- pi = 1 - pi ;
140+ pi = ( I2 > 0 ) ? 1 - pi : 0 ;
141141 event_rate += (* rate = alpha * pi ); rate ++ ;
142142 * logpi = log (pi )- log (ell2 ); logpi ++ ;
143143 assert (R_FINITE (event_rate ));
0 commit comments