]> henry.ined.fr Git - .git/commitdiff
*** empty log message ***
authorAgnès Lièvre <agnes.lievre@education.gouv.fr>
Fri, 5 Apr 2002 15:45:00 +0000 (15:45 +0000)
committerAgnès Lièvre <agnes.lievre@education.gouv.fr>
Fri, 5 Apr 2002 15:45:00 +0000 (15:45 +0000)
src/imach.c

index 28c1ec7335fa959c02f284e0a42fbf6750b27ac2..5ba0474b40bcacc6af8a5f747863bce0efbe8940 100644 (file)
@@ -14,7 +14,7 @@
   model. More health states you consider, more time is necessary to reach the\r
   Maximum Likelihood of the parameters involved in the model.  The\r
   simplest model is the multinomial logistic model where pij is the\r
-  probabibility to be observed in state j at the second wave\r
+  probability to be observed in state j at the second wave\r
   conditional to be observed in state i at the first wave. Therefore\r
   the model is: log(pij/pii)= aij + bij*age+ cij*sex + etc , where\r
   'age' is age and 'sex' is a covariate. If you want to have a more\r
@@ -1503,7 +1503,7 @@ void tricode(int *Tvar, int **nbcode, int imx)
       for (k=0; k<=19; k++) {\r
        if (Ndum[k] != 0) {\r
          nbcode[Tvar[j]][ij]=k; \r
-         /*     printf("nbcodeaaaaaaaaaaa=%d Tvar[j]=%d ij=%d j=%d",nbcode[Tvar[j]][ij],Tvar[j],ij,j);*/\r
+         \r
          ij++;\r
        }\r
        if (ij > ncodemax[j]) break; \r
@@ -1585,6 +1585,9 @@ void evsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, in
     /* Computed by stepm unit matrices, product of hstepm matrices, stored\r
        in an array of nhstepm length: nhstepm=10, hstepm=4, stepm=6 months */\r
     hpxij(p3mat,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, ij);  \r
+  \r
+    /*for (h=0, eij[i][j][(int)age]=0; h<=nhstepm-1; h++) printf("%f %.5f\n", age*12+h, p3mat[1][1][h]);*/\r
+\r
     hf=hstepm*stepm/YEARM;  /* Duration of hstepm expressed in year unit. */\r
     for(i=1; i<=nlstate;i++)\r
       for(j=1; j<=nlstate;j++)\r
@@ -1822,9 +1825,9 @@ void varprevlim(char fileres[], double **varpl, double **matcov, double x[], dou
 }\r
 \r
 /************ Variance of one-step probabilities  ******************/\r
-void varprob(char fileres[], double **matcov, double x[], double delti[], int nlstate, double bage, double fage, int ij)\r
+void varprob(char fileres[], double **matcov, double x[], double delti[], int nlstate, double bage, double fage, int ij, int *Tvar, int **nbcode, int *ncodemax)\r
 {\r
-  int i, j;\r
+  int i, j, i1, k1, j1, z1;\r
   int k=0, cptcode;\r
   double **dnewm,**doldm;\r
   double *xp;\r
@@ -1847,83 +1850,101 @@ void varprob(char fileres[], double **matcov, double x[], double delti[], int nl
   doldm=matrix(1,(nlstate+ndeath)*(nlstate+ndeath),1,(nlstate+ndeath)*(nlstate+ndeath));\r
   \r
   cov[1]=1;\r
-  for (age=bage; age<=fage; age ++){ \r
-    cov[2]=age;\r
-    gradg=matrix(1,npar,1,9);\r
-    trgradg=matrix(1,9,1,npar);\r
-    gp=vector(1,(nlstate+ndeath)*(nlstate+ndeath));\r
-    gm=vector(1,(nlstate+ndeath)*(nlstate+ndeath));\r
+  j=cptcoveff;\r
+  if (cptcovn<1) {j=1;ncodemax[1]=1;}\r
+  j1=0;\r
+  for(k1=1; k1<=1;k1++){\r
+    for(i1=1; i1<=ncodemax[k1];i1++){ \r
+    j1++;\r
+\r
+    if  (cptcovn>0) {\r
+      fprintf(ficresprob, "\n#********** Variable "); \r
+      for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprob, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtab[j1][z1]]);\r
+      fprintf(ficresprob, "**********\n#");\r
+    }\r
     \r
-    for(theta=1; theta <=npar; theta++){\r
-      for(i=1; i<=npar; i++)\r
-       xp[i] = x[i] + (i==theta ?delti[theta]:0);\r
-     \r
-      pmij(pmmij,cov,ncovmodel,xp,nlstate);\r
-   \r
-      k=0;\r
-      for(i=1; i<= (nlstate+ndeath); i++){\r
-       for(j=1; j<=(nlstate+ndeath);j++){\r
-          k=k+1;\r
-         gp[k]=pmmij[i][j];\r
+      for (age=bage; age<=fage; age ++){ \r
+       cov[2]=age;\r
+       for (k=1; k<=cptcovn;k++) {\r
+         cov[2+k]=nbcode[Tvar[k]][codtab[j1][Tvar[k]]];\r
+         \r
        }\r
-      }\r
-\r
-      for(i=1; i<=npar; i++)\r
-       xp[i] = x[i] - (i==theta ?delti[theta]:0);\r
+       for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2];\r
+       for (k=1; k<=cptcovprod;k++)\r
+         cov[2+Tprod[k]]=nbcode[Tvard[k][1]][codtab[ij][Tvard[k][1]]]*nbcode[Tvard[k][2]][codtab[ij][Tvard[k][2]]];\r
+       \r
+       gradg=matrix(1,npar,1,9);\r
+       trgradg=matrix(1,9,1,npar);\r
+       gp=vector(1,(nlstate+ndeath)*(nlstate+ndeath));\r
+       gm=vector(1,(nlstate+ndeath)*(nlstate+ndeath));\r
     \r
-\r
-      pmij(pmmij,cov,ncovmodel,xp,nlstate);\r
-      k=0;\r
-      for(i=1; i<=(nlstate+ndeath); i++){\r
-       for(j=1; j<=(nlstate+ndeath);j++){\r
-         k=k+1;\r
-         gm[k]=pmmij[i][j];\r
-       }\r
-      }\r
+       for(theta=1; theta <=npar; theta++){\r
+         for(i=1; i<=npar; i++)\r
+           xp[i] = x[i] + (i==theta ?delti[theta]:0);\r
+         \r
+         pmij(pmmij,cov,ncovmodel,xp,nlstate);\r
+         \r
+         k=0;\r
+         for(i=1; i<= (nlstate+ndeath); i++){\r
+           for(j=1; j<=(nlstate+ndeath);j++){\r
+             k=k+1;\r
+             gp[k]=pmmij[i][j];\r
+           }\r
+         }\r
+         \r
+         for(i=1; i<=npar; i++)\r
+           xp[i] = x[i] - (i==theta ?delti[theta]:0);\r
+    \r
+         pmij(pmmij,cov,ncovmodel,xp,nlstate);\r
+         k=0;\r
+         for(i=1; i<=(nlstate+ndeath); i++){\r
+           for(j=1; j<=(nlstate+ndeath);j++){\r
+             k=k+1;\r
+             gm[k]=pmmij[i][j];\r
+           }\r
+         }\r
      \r
-       for(i=1; i<= (nlstate+ndeath)*(nlstate+ndeath); i++) \r
-          gradg[theta][i]=(gp[i]-gm[i])/2./delti[theta];  \r
-    }\r
-\r
-     for(j=1; j<=(nlstate+ndeath)*(nlstate+ndeath);j++)\r
-      for(theta=1; theta <=npar; theta++)\r
-      trgradg[j][theta]=gradg[theta][j];\r
\r
-     matprod2(dnewm,trgradg,1,9,1,npar,1,npar,matcov);\r
-     matprod2(doldm,dnewm,1,9,1,npar,1,9,gradg);\r
-\r
-     pmij(pmmij,cov,ncovmodel,x,nlstate);\r
+         for(i=1; i<= (nlstate+ndeath)*(nlstate+ndeath); i++) \r
+           gradg[theta][i]=(gp[i]-gm[i])/2./delti[theta];  \r
+       }\r
 \r
-     k=0;\r
-     for(i=1; i<=(nlstate+ndeath); i++){\r
-       for(j=1; j<=(nlstate+ndeath);j++){\r
-        k=k+1;\r
-        gm[k]=pmmij[i][j];\r
+       for(j=1; j<=(nlstate+ndeath)*(nlstate+ndeath);j++)\r
+         for(theta=1; theta <=npar; theta++)\r
+           trgradg[j][theta]=gradg[theta][j];\r
+       \r
+       matprod2(dnewm,trgradg,1,9,1,npar,1,npar,matcov);\r
+       matprod2(doldm,dnewm,1,9,1,npar,1,9,gradg);\r
+       \r
+       pmij(pmmij,cov,ncovmodel,x,nlstate);\r
+       \r
+       k=0;\r
+       for(i=1; i<=(nlstate+ndeath); i++){\r
+         for(j=1; j<=(nlstate+ndeath);j++){\r
+           k=k+1;\r
+           gm[k]=pmmij[i][j];\r
+         }\r
        }\r
-     }\r
      \r
      /*printf("\n%d ",(int)age);\r
      for (i=1; i<=(nlstate+ndeath)*(nlstate+ndeath-1);i++){\r
-       \r
-\r
        printf("%e [%e ;%e] ",gm[i],gm[i]-2*sqrt(doldm[i][i]),gm[i]+2*sqrt(doldm[i][i]));\r
      }*/\r
 \r
-  fprintf(ficresprob,"\n%d ",(int)age);\r
-\r
-  for (i=1; i<=(nlstate+ndeath)*(nlstate+ndeath-1);i++){\r
-    if (i== 2) fprintf(ficresprob,"%.3e %.3e ",gm[i],doldm[i][i]);\r
-if (i== 4) fprintf(ficresprob,"%.3e %.3e ",gm[i],doldm[i][i]);\r
-  }\r
+       fprintf(ficresprob,"\n%d ",(int)age);\r
 \r
+       for (i=1; i<=(nlstate+ndeath)*(nlstate+ndeath-1);i++)\r
+         fprintf(ficresprob,"%.3e (%.3e) ",gm[i],doldm[i][i]);\r
+  \r
+      }\r
+    }\r
     free_vector(gp,1,(nlstate+ndeath)*(nlstate+ndeath));\r
     free_vector(gm,1,(nlstate+ndeath)*(nlstate+ndeath));\r
     free_matrix(trgradg,1,(nlstate+ndeath)*(nlstate+ndeath),1,npar);\r
     free_matrix(gradg,1,(nlstate+ndeath)*(nlstate+ndeath),1,npar);\r
-}\r
- free_vector(xp,1,npar);\r
-fclose(ficresprob);\r
-\r
+  }\r
 free_vector(xp,1,npar);\r
+  fclose(ficresprob);\r
+  \r
 }\r
 \r
 /******************* Printing html file ***********/\r
@@ -2751,11 +2772,10 @@ while((c=getc(ficpar))=='#' && c!= EOF){
     if ((s[2][i]==3) && (s[3][i]==2)) s[3][i]=3;\r
     if ((s[3][i]==3) && (s[4][i]==2)) s[4][i]=3;\r
     }*/\r
\r
-  /* for (i=1; i<=imx; i++){\r
+   /*  for (i=1; i<=imx; i++){\r
      if (s[4][i]==9)  s[4][i]=-1; \r
-     printf("%d %.lf %.lf %.lf %.lf/%.lf %.lf/%.lf %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d\n",num[i],(covar[1][i]), (covar[2][i]), (weight[i]), (moisnais[i]), (annais[i]), (moisdc[i]), (andc[i]), (mint[1][i]), (anint[1][i]), (s[1][i]),  (mint[2][i]), (anint[2][i]), (s[2][i]),  (mint[3][i]), (anint[3][i]), (s[3][i]),  (mint[4][i]), (anint[4][i]), (s[4][i]));}\r
-  */\r
+     printf("%d %.lf %.lf %.lf %.lf/%.lf %.lf/%.lf %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d\n",num[i],(covar[1][i]), (covar[2][i]), (weight[i]), (moisnais[i]), (annais[i]), (moisdc[i]), (andc[i]), (mint[1][i]), (anint[1][i]), (s[1][i]),  (mint[2][i]), (anint[2][i]), (s[2][i]),  (mint[3][i]), (anint[3][i]), (s[3][i]),  (mint[4][i]), (anint[4][i]), (s[4][i]));}*/\r
+  \r
  \r
   /* Calculation of the number of parameter from char model*/\r
   Tvar=ivector(1,15); \r
@@ -3220,7 +3240,7 @@ while((c=getc(ficpar))=='#' && c!= EOF){
     }\r
   }\r
 \r
-  /* varprob(fileres, matcov, p, delti, nlstate, (int) bage, (int) fage,k);*/\r
+  varprob(fileres, matcov, p, delti, nlstate, (int) bage, (int) fage,k,Tvar,nbcode, ncodemax);\r
 \r
   fclose(ficrespij);\r
 \r