We decided to release, in June 2003, a version 0.8a1 (May 2002) of

former version 0.8a to correct a bug highlighted by Chris Jackson.
It can be useful until version 0.9x will be released.

In version 0.8a, likelihood was computed as if date of death was unknown. Death
was treated as any other health state, i.e the date of the interview
describes the actual state and not the date of a change in health
state.

It means that in the early version of IMaCh we were considering that
at each interview the current state was recorded (healthy, disable or
death).  And IMaCh was correct. But, when we introduced the exact date
of death then we should have modified the contribution of an exact
(more precise) date of death to the likelihood. This new contribution
is smaller and very dependent of the step unit stepm. It is no more
the probability to die between last interview and month of death
(which is the way version 0.8a of IMaCh was doing) but the probability
to survive from last interview up to one month before death multiplied
by the probability to die within a month. Thanks to Chris Jackson for
correcting this bug.  Former versions increased mortality
artificially. The bad side is that we add another loop which slows
down the processing. The difference can be up to 10% lower mortality.

diff --git a/src/imach.c b/src/imach.c
index 636f597133a74539bb4ad4d78b48b0da79c43281..3efca52db84fd81501922a93bd1e7d7d6c4dd74f 100644 (file)
--- a/src/imach.c
+++ b/src/imach.c
@@ -870,6 +870,7 @@ double func( double *x)
   double **out;\r
   double sw; /* Sum of weights */\r
   double lli; /* Individual log likelihood */\r
+  int s1, s2;\r
   long ipmx;\r
   /*extern weight */\r
   /* We are differentiating ll according to initial status */\r
@@ -884,7 +885,10 @@ double func( double *x)
     for (k=1; k<=cptcovn;k++) cov[2+k]=covar[Tvar[k]][i];\r
     for(mi=1; mi<= wav[i]-1; mi++){\r
       for (ii=1;ii<=nlstate+ndeath;ii++)\r
-       for (j=1;j<=nlstate+ndeath;j++) oldm[ii][j]=(ii==j ? 1.0 : 0.0);\r
+       for (j=1;j<=nlstate+ndeath;j++){\r
+         oldm[ii][j]=(ii==j ? 1.0 : 0.0);\r
+         savm[ii][j]=(ii==j ? 1.0 : 0.0);\r
+       }\r
       for(d=0; d<dh[mi][i]; d++){\r
        newm=savm;\r
        cov[2]=agev[mw[mi][i]][i]+d*stepm/YEARM;\r
@@ -900,17 +904,47 @@ double func( double *x)
        \r
       } /* end mult */\r
       \r
-      lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);\r
-      /* printf(" %f ",out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/\r
+      s1=s[mw[mi][i]][i];\r
+      s2=s[mw[mi+1][i]][i];\r
+      if( s2 > nlstate){ \r
+       /* i.e. if s2 is a death state and if the date of death is known then the contribution\r
+          to the likelihood is the probability to die between last step unit time and current \r
+          step unit time, which is also the differences between probability to die before dh \r
+          and probability to die before dh-stepm . \r
+          In version up to 0.92 likelihood was computed\r
+          as if date of death was unknown. Death was treated as any other\r
+          health state: the date of the interview describes the actual state\r
+          and not the date of a change in health state. The former idea was\r
+          to consider that at each interview the state was recorded\r
+          (healthy, disable or death) and IMaCh was corrected; but when we\r
+          introduced the exact date of death then we should have modified\r
+          the contribution of an exact death to the likelihood. This new\r
+          contribution is smaller and very dependent of the step unit\r
+          stepm. It is no more the probability to die between last interview\r
+          and month of death but the probability to survive from last\r
+          interview up to one month before death multiplied by the\r
+          probability to die within a month. Thanks to Chris\r
+          Jackson for correcting this bug.  Former versions increased\r
+          mortality artificially. The bad side is that we add another loop\r
+          which slows down the processing. The difference can be up to 10%\r
+          lower mortality.\r
+       */\r
+       lli=log(out[s1][s2] - savm[s1][s2]);\r
+      }else{\r
+       lli=log(out[s1][s2]); /* or     lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]); */\r
+       /* printf(" %f ",out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/\r
+      }\r
       ipmx +=1;\r
       sw += weight[i];\r
       ll[s[mw[mi][i]][i]] += 2*weight[i]*lli;\r
+      /*printf("i=%6d s1=%1d s2=%1d mi=%1d mw=%1d lli=%10.6f w=%6.4f out=%10.6f sav=%10.6f\n",i,s1,s2,mi,mw[mi][i],lli,weight[i],out[s1][s2],savm[s1][s2]);*/\r
     } /* end of wave */\r
   } /* end of individual */\r
 \r
   for(k=1,l=0.; k<=nlstate; k++) l += ll[k];\r
   /* printf("l1=%f l2=%f ",ll[1],ll[2]); */\r
   l= l*ipmx/sw; /* To get the same order of magnitude as if weight=1 for every body */\r
+  /*exit(0);*/\r
   return -l;\r
 }\r
 \r
@@ -2645,7 +2679,7 @@ int main(int argc, char *argv[])
   double dateprev1, dateprev2,jproj1,mproj1,anproj1,jproj2,mproj2,anproj2;\r
   \r
 \r
-  char version[80]="Imach version 0.8a, May 2002, INED-EUROREVES ";\r
+  char version[80]="Imach version 0.8a1, June 2003, INED-EUROREVES ";\r
   char *alph[]={"a","a","b","c","d","e"}, str[4];\r
 \r
 \r