|
version 1.221, 2016/02/15 23:35:36
|
version 1.222, 2016/02/17 08:14:50
|
|
Line 1
|
Line 1
|
| /* $Id$ |
/* $Id$ |
| $State$ |
$State$ |
| $Log$ |
$Log$ |
| |
Revision 1.222 2016/02/17 08:14:50 brouard |
| |
Summary: Probably last 0.98 stable version 0.98r6 |
| |
|
| Revision 1.221 2016/02/15 23:35:36 brouard |
Revision 1.221 2016/02/15 23:35:36 brouard |
| Summary: minor bug |
Summary: minor bug |
| |
|
|
Line 2398 double **pmij(double **ps, double *cov,
|
Line 2401 double **pmij(double **ps, double *cov,
|
| /* double **bmij(double **ps, double *cov, int ncovmodel, double *x, int nlstate, double ***prevacurrent, double ***dnewm, double **doldm, double **dsavm, int ij ) */ |
/* double **bmij(double **ps, double *cov, int ncovmodel, double *x, int nlstate, double ***prevacurrent, double ***dnewm, double **doldm, double **dsavm, int ij ) */ |
| double **bmij(double **ps, double *cov, int ncovmodel, double *x, int nlstate, double ***prevacurrent, int ij ) |
double **bmij(double **ps, double *cov, int ncovmodel, double *x, int nlstate, double ***prevacurrent, int ij ) |
| { |
{ |
| /* Computes the backward probability at age agefin and covariate ij |
/* Computes the backward probability at age agefin and covariate ij |
| * and returns in **ps as well as **bmij. |
* and returns in **ps as well as **bmij. |
| */ |
*/ |
| int i, ii, j,k; |
int i, ii, j,k; |
| |
|
| double **out, **pmij(); |
double **out, **pmij(); |
| double sumnew=0.; |
double sumnew=0.; |
| double agefin; |
double agefin; |
| |
|
| double **dnewm, **dsavm, **doldm; |
double **dnewm, **dsavm, **doldm; |
| double **bbmij; |
double **bbmij; |
| |
|
| doldm=ddoldms; /* global pointers */ |
doldm=ddoldms; /* global pointers */ |
| dnewm=ddnewms; |
dnewm=ddnewms; |
| dsavm=ddsavms; |
dsavm=ddsavms; |
| |
|
| agefin=cov[2]; |
agefin=cov[2]; |
| /* bmij *//* age is cov[2], ij is included in cov, but we need for |
/* bmij *//* age is cov[2], ij is included in cov, but we need for |
| the observed prevalence (with this covariate ij) */ |
the observed prevalence (with this covariate ij) */ |
| dsavm=pmij(pmmij,cov,ncovmodel,x,nlstate); |
dsavm=pmij(pmmij,cov,ncovmodel,x,nlstate); |
| /* We do have the matrix Px in savm and we need pij */ |
/* We do have the matrix Px in savm and we need pij */ |
| for (j=1;j<=nlstate+ndeath;j++){ |
for (j=1;j<=nlstate+ndeath;j++){ |
| sumnew=0.; /* w1 p11 + w2 p21 only on live states */ |
sumnew=0.; /* w1 p11 + w2 p21 only on live states */ |
| for (ii=1;ii<=nlstate;ii++){ |
for (ii=1;ii<=nlstate;ii++){ |
| sumnew+=dsavm[ii][j]*prevacurrent[(int)agefin][ii][ij]; |
sumnew+=dsavm[ii][j]*prevacurrent[(int)agefin][ii][ij]; |
| } /* sumnew is (N11+N21)/N..= N.1/N.. = sum on i of w_i pij */ |
} /* sumnew is (N11+N21)/N..= N.1/N.. = sum on i of w_i pij */ |
| for (ii=1;ii<=nlstate+ndeath;ii++){ |
for (ii=1;ii<=nlstate+ndeath;ii++){ |
| if(sumnew >= 1.e-10){ |
if(sumnew >= 1.e-10){ |
| /* if(agefin >= agemaxpar && agefin <= agemaxpar+stepm/YEARM){ */ |
/* if(agefin >= agemaxpar && agefin <= agemaxpar+stepm/YEARM){ */ |
| /* doldm[ii][j]=(ii==j ? 1./sumnew : 0.0); */ |
/* doldm[ii][j]=(ii==j ? 1./sumnew : 0.0); */ |
| /* }else if(agefin >= agemaxpar+stepm/YEARM){ */ |
/* }else if(agefin >= agemaxpar+stepm/YEARM){ */ |
| /* doldm[ii][j]=(ii==j ? 1./sumnew : 0.0); */ |
/* doldm[ii][j]=(ii==j ? 1./sumnew : 0.0); */ |
| /* }else */ |
/* }else */ |
| doldm[ii][j]=(ii==j ? 1./sumnew : 0.0); |
doldm[ii][j]=(ii==j ? 1./sumnew : 0.0); |
| }else{ |
}else{ |
| printf("ii=%d, i=%d, doldm=%lf dsavm=%lf, probs=%lf, sumnew=%lf,agefin=%d\n",ii,j,doldm[ii][j],dsavm[ii][j],prevacurrent[(int)agefin][ii][ij],sumnew, (int)agefin); |
printf("ii=%d, i=%d, doldm=%lf dsavm=%lf, probs=%lf, sumnew=%lf,agefin=%d\n",ii,j,doldm[ii][j],dsavm[ii][j],prevacurrent[(int)agefin][ii][ij],sumnew, (int)agefin); |
| } |
} |
| } /*End ii */ |
} /*End ii */ |
| } /* End j, At the end doldm is diag[1/(w_1p1i+w_2 p2i)] */ |
} /* End j, At the end doldm is diag[1/(w_1p1i+w_2 p2i)] */ |
| /* left Product of this diag matrix by dsavm=Px (newm=dsavm*doldm) */ |
/* left Product of this diag matrix by dsavm=Px (newm=dsavm*doldm) */ |
| bbmij=matprod2(dnewm, dsavm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, doldm); /* Bug Valgrind */ |
bbmij=matprod2(dnewm, dsavm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, doldm); /* Bug Valgrind */ |
| /* dsavm=doldm; /\* dsavm is now diag [1/(w_1p1i+w_2 p2i)] but can be overwritten*\/ */ |
/* dsavm=doldm; /\* dsavm is now diag [1/(w_1p1i+w_2 p2i)] but can be overwritten*\/ */ |
| /* doldm=dnewm; /\* doldm is now Px * diag [1/(w_1p1i+w_2 p2i)] *\/ */ |
/* doldm=dnewm; /\* doldm is now Px * diag [1/(w_1p1i+w_2 p2i)] *\/ */ |
| /* dnewm=dsavm; /\* doldm is now Px * diag [1/(w_1p1i+w_2 p2i)] *\/ */ |
/* dnewm=dsavm; /\* doldm is now Px * diag [1/(w_1p1i+w_2 p2i)] *\/ */ |
| /* left Product of this matrix by diag matrix of prevalences (savm) */ |
/* left Product of this matrix by diag matrix of prevalences (savm) */ |
| for (j=1;j<=nlstate+ndeath;j++){ |
for (j=1;j<=nlstate+ndeath;j++){ |
| for (ii=1;ii<=nlstate+ndeath;ii++){ |
for (ii=1;ii<=nlstate+ndeath;ii++){ |
| dsavm[ii][j]=(ii==j ? prevacurrent[(int)agefin][ii][ij] : 0.0); |
dsavm[ii][j]=(ii==j ? prevacurrent[(int)agefin][ii][ij] : 0.0); |
| } |
} |
| } /* End j, At the end oldm is diag[1/(w_1p1i+w_2 p2i)] */ |
} /* End j, At the end oldm is diag[1/(w_1p1i+w_2 p2i)] */ |
| ps=matprod2(doldm, dsavm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, dnewm); /* Bug Valgrind */ |
ps=matprod2(doldm, dsavm,1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, dnewm); /* Bug Valgrind */ |
| /* newm or out is now diag[w_i] * Px * diag [1/(w_1p1i+w_2 p2i)] */ |
/* newm or out is now diag[w_i] * Px * diag [1/(w_1p1i+w_2 p2i)] */ |
| /* end bmij */ |
/* end bmij */ |
| return ps; |
return ps; |
| } |
} |
| /*************** transition probabilities ***************/ |
/*************** transition probabilities ***************/ |
| |
|
|
Line 2657 double ***hpxij(double ***po, int nhstep
|
Line 2660 double ***hpxij(double ***po, int nhstep
|
| |
|
| /************* Higher Back Matrix Product ***************/ |
/************* Higher Back Matrix Product ***************/ |
| /* double ***hbxij(double ***po, int nhstepm, double age, int hstepm, double *x, double ***prevacurrent, int nlstate, int stepm, double **oldm, double **savm, double **dnewm, double **doldm, double **dsavm, int ij ) */ |
/* double ***hbxij(double ***po, int nhstepm, double age, int hstepm, double *x, double ***prevacurrent, int nlstate, int stepm, double **oldm, double **savm, double **dnewm, double **doldm, double **dsavm, int ij ) */ |
| double ***hbxij(double ***po, int nhstepm, double age, int hstepm, double *x, double ***prevacurrent, int nlstate, int stepm, int ij ) |
double ***hbxij(double ***po, int nhstepm, double age, int hstepm, double *x, double ***prevacurrent, int nlstate, int stepm, int ij ) |
| { |
{ |
| /* Computes the transition matrix starting at age 'age' over |
/* Computes the transition matrix starting at age 'age' over |
| 'nhstepm*hstepm*stepm' months (i.e. until |
'nhstepm*hstepm*stepm' months (i.e. until |
|
Line 2669 double ***hpxij(double ***po, int nhstep
|
Line 2672 double ***hpxij(double ***po, int nhstep
|
| Model is determined by parameters x and covariates have to be |
Model is determined by parameters x and covariates have to be |
| included manually here. |
included manually here. |
| |
|
| */ |
*/ |
| |
|
| int i, j, d, h, k; |
int i, j, d, h, k; |
| double **out, cov[NCOVMAX+1]; |
double **out, cov[NCOVMAX+1]; |
| double **newm; |
double **newm; |
| double agexact; |
double agexact; |
| double agebegin, ageend; |
double agebegin, ageend; |
| double **oldm, **savm; |
double **oldm, **savm; |
| |
|
| oldm=oldms;savm=savms; |
oldm=oldms;savm=savms; |
| /* Hstepm could be zero and should return the unit matrix */ |
/* Hstepm could be zero and should return the unit matrix */ |
| for (i=1;i<=nlstate+ndeath;i++) |
for (i=1;i<=nlstate+ndeath;i++) |
| for (j=1;j<=nlstate+ndeath;j++){ |
for (j=1;j<=nlstate+ndeath;j++){ |
|
Line 2695 double ***hpxij(double ***po, int nhstep
|
Line 2698 double ***hpxij(double ***po, int nhstep
|
| /* agexact=age+((h-1)*hstepm + (d-1))*stepm/YEARM; /\* age just before transition *\/ */ |
/* agexact=age+((h-1)*hstepm + (d-1))*stepm/YEARM; /\* age just before transition *\/ */ |
| cov[2]=agexact; |
cov[2]=agexact; |
| if(nagesqr==1) |
if(nagesqr==1) |
| cov[3]= agexact*agexact; |
cov[3]= agexact*agexact; |
| for (k=1; k<=cptcovn;k++) |
for (k=1; k<=cptcovn;k++) |
| cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,k)]; |
cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,k)]; |
| /* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,Tvar[k])]; */ |
/* cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(ij,Tvar[k])]; */ |
| for (k=1; k<=cptcovage;k++) /* Should start at cptcovn+1 */ |
for (k=1; k<=cptcovage;k++) /* Should start at cptcovn+1 */ |
| /* cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */ |
/* cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */ |
| cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; |
cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; |
| /* cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,Tvar[Tage[k]])]*cov[2]; */ |
/* cov[2+nagesqr+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,Tvar[Tage[k]])]*cov[2]; */ |
| for (k=1; k<=cptcovprod;k++) /* Useless because included in cptcovn */ |
for (k=1; k<=cptcovprod;k++) /* Useless because included in cptcovn */ |
| cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)]; |
cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)]; |
| /* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,Tvard[k][1])]*nbcode[Tvard[k][2]][codtabm(ij,Tvard[k][2])]; */ |
/* cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,Tvard[k][1])]*nbcode[Tvard[k][2]][codtabm(ij,Tvard[k][2])]; */ |
| |
|
| |
|
| /*printf("hxi cptcov=%d cptcode=%d\n",cptcov,cptcode);*/ |
/*printf("hxi cptcov=%d cptcode=%d\n",cptcov,cptcode);*/ |
| /*printf("h=%d d=%d age=%f cov=%f\n",h,d,age,cov[2]);*/ |
/*printf("h=%d d=%d age=%f cov=%f\n",h,d,age,cov[2]);*/ |
| /* Careful transposed matrix */ |
/* Careful transposed matrix */ |
| /* age is in cov[2] */ |
/* age is in cov[2] */ |
| /* out=matprod2(newm, bmij(pmmij,cov,ncovmodel,x,nlstate,prevacurrent, dnewm, doldm, dsavm,ij),\ */ |
/* out=matprod2(newm, bmij(pmmij,cov,ncovmodel,x,nlstate,prevacurrent, dnewm, doldm, dsavm,ij),\ */ |
| /* 1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, oldm); */ |
/* 1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, oldm); */ |
| out=matprod2(newm, bmij(pmmij,cov,ncovmodel,x,nlstate,prevacurrent,ij),\ |
out=matprod2(newm, bmij(pmmij,cov,ncovmodel,x,nlstate,prevacurrent,ij),\ |
| 1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, oldm); |
1,nlstate+ndeath,1,nlstate+ndeath,1,nlstate+ndeath, oldm); |
| /* if((int)age == 70){ */ |
/* if((int)age == 70){ */ |
| /* printf(" Backward hbxij age=%d agexact=%f d=%d nhstepm=%d hstepm=%d\n", (int) age, agexact, d, nhstepm, hstepm); */ |
/* printf(" Backward hbxij age=%d agexact=%f d=%d nhstepm=%d hstepm=%d\n", (int) age, agexact, d, nhstepm, hstepm); */ |
| /* for(i=1; i<=nlstate+ndeath; i++) { */ |
/* for(i=1; i<=nlstate+ndeath; i++) { */ |
|
Line 2735 double ***hpxij(double ***po, int nhstep
|
Line 2738 double ***hpxij(double ***po, int nhstep
|
| } |
} |
| for(i=1; i<=nlstate+ndeath; i++) |
for(i=1; i<=nlstate+ndeath; i++) |
| for(j=1;j<=nlstate+ndeath;j++) { |
for(j=1;j<=nlstate+ndeath;j++) { |
| po[i][j][h]=newm[i][j]; |
po[i][j][h]=newm[i][j]; |
| /*if(h==nhstepm) printf("po[%d][%d][%d]=%f ",i,j,h,po[i][j][h]);*/ |
/*if(h==nhstepm) printf("po[%d][%d][%d]=%f ",i,j,h,po[i][j][h]);*/ |
| } |
} |
| /*printf("h=%d ",h);*/ |
/*printf("h=%d ",h);*/ |
| } /* end h */ |
} /* end h */ |
| /* printf("\n H=%d \n",h); */ |
/* printf("\n H=%d \n",h); */ |
| return po; |
return po; |
| } |
} |
| |
|
|
Line 3995 Title=%s <br>Datafile=%s Firstpass=%d La
|
Line 3998 Title=%s <br>Datafile=%s Firstpass=%d La
|
| } |
} |
| |
|
| /************ Prevalence ********************/ |
/************ Prevalence ********************/ |
| void prevalence(double ***probs, double agemin, double agemax, int **s, double **agev, int nlstate, int imx, int *Tvar, int **nbcode, int *ncodemax,double **mint,double **anint, double dateprev1,double dateprev2, int firstpass, int lastpass) |
void prevalence(double ***probs, double agemin, double agemax, int **s, double **agev, int nlstate, int imx, int *Tvar, int **nbcode, int *ncodemax,double **mint,double **anint, double dateprev1,double dateprev2, int firstpass, int lastpass) |
| { |
{ |
| /* Compute observed prevalence between dateprev1 and dateprev2 by counting the number of people |
/* Compute observed prevalence between dateprev1 and dateprev2 by counting the number of people |
| in each health status at the date of interview (if between dateprev1 and dateprev2). |
in each health status at the date of interview (if between dateprev1 and dateprev2). |
| We still use firstpass and lastpass as another selection. |
We still use firstpass and lastpass as another selection. |
| */ |
*/ |
| |
|
| int i, m, jk, j1, bool, z1,j; |
int i, m, jk, j1, bool, z1,j; |
| int mi; /* Effective wave */ |
int mi; /* Effective wave */ |
| int iage; |
int iage; |
| double agebegin, ageend; |
double agebegin, ageend; |
| |
|
| double **prop; |
double **prop; |
| double posprop; |
double posprop; |
| double y2; /* in fractional years */ |
double y2; /* in fractional years */ |
| int iagemin, iagemax; |
int iagemin, iagemax; |
| int first; /** to stop verbosity which is redirected to log file */ |
int first; /** to stop verbosity which is redirected to log file */ |
| |
|
| iagemin= (int) agemin; |
iagemin= (int) agemin; |
| iagemax= (int) agemax; |
iagemax= (int) agemax; |
| /*pp=vector(1,nlstate);*/ |
/*pp=vector(1,nlstate);*/ |
| prop=matrix(1,nlstate,iagemin-AGEMARGE,iagemax+3+AGEMARGE); |
prop=matrix(1,nlstate,iagemin-AGEMARGE,iagemax+3+AGEMARGE); |
| /* freq=ma3x(-1,nlstate+ndeath,-1,nlstate+ndeath,iagemin,iagemax+3);*/ |
/* freq=ma3x(-1,nlstate+ndeath,-1,nlstate+ndeath,iagemin,iagemax+3);*/ |
| j1=0; |
j1=0; |
| |
|
| /*j=cptcoveff;*/ |
/*j=cptcoveff;*/ |
| if (cptcovn<1) {j=1;ncodemax[1]=1;} |
if (cptcovn<1) {j=1;ncodemax[1]=1;} |
| |
|
| first=1; |
first=1; |
| for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of covariate */ |
for(j1=1; j1<= (int) pow(2,cptcoveff);j1++){ /* For each combination of covariate */ |
| for (i=1; i<=nlstate; i++) |
for (i=1; i<=nlstate; i++) |
| for(iage=iagemin-AGEMARGE; iage <= iagemax+3+AGEMARGE; iage++) |
for(iage=iagemin-AGEMARGE; iage <= iagemax+3+AGEMARGE; iage++) |
| prop[i][iage]=0.0; |
prop[i][iage]=0.0; |
| |
|
| for (i=1; i<=imx; i++) { /* Each individual */ |
for (i=1; i<=imx; i++) { /* Each individual */ |
| bool=1; |
bool=1; |
| if (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */ |
if (cptcovn>0) { /* Filter is here: Must be looked at for model=V1+V2+V3+V4 */ |
| for (z1=1; z1<=cptcoveff; z1++) /* For each covariate, look at the value for individual i and checks if it is equal to the corresponding value of this covariate according to current combination j1*/ |
for (z1=1; z1<=cptcoveff; z1++) /* For each covariate, look at the value for individual i and checks if it is equal to the corresponding value of this covariate according to current combination j1*/ |
| if (covar[Tvaraff[z1]][i]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]) |
if (covar[Tvaraff[z1]][i]!= nbcode[Tvaraff[z1]][codtabm(j1,z1)]) |
| bool=0; |
bool=0; |
| } |
} |
| if (bool==1) { /* For this combination of covariates values, this individual fits */ |
if (bool==1) { /* For this combination of covariates values, this individual fits */ |
| /* for(m=firstpass; m<=lastpass; m++){/\* Other selection (we can limit to certain interviews*\/ */ |
/* for(m=firstpass; m<=lastpass; m++){/\* Other selection (we can limit to certain interviews*\/ */ |
| for(mi=1; mi<wav[i];mi++){ |
for(mi=1; mi<wav[i];mi++){ |
| m=mw[mi][i]; |
m=mw[mi][i]; |
| agebegin=agev[m][i]; /* Age at beginning of wave before transition*/ |
agebegin=agev[m][i]; /* Age at beginning of wave before transition*/ |
| /* ageend=agev[m][i]+(dh[m][i])*stepm/YEARM; /\* Age at end of wave and transition *\/ */ |
/* ageend=agev[m][i]+(dh[m][i])*stepm/YEARM; /\* Age at end of wave and transition *\/ */ |
| if(m >=firstpass && m <=lastpass){ |
if(m >=firstpass && m <=lastpass){ |
| y2=anint[m][i]+(mint[m][i]/12.); /* Fractional date in year */ |
y2=anint[m][i]+(mint[m][i]/12.); /* Fractional date in year */ |
| if ((y2>=dateprev1) && (y2<=dateprev2)) { /* Here is the main selection (fractional years) */ |
if ((y2>=dateprev1) && (y2<=dateprev2)) { /* Here is the main selection (fractional years) */ |
| if(agev[m][i]==0) agev[m][i]=iagemax+1; |
if(agev[m][i]==0) agev[m][i]=iagemax+1; |
| if(agev[m][i]==1) agev[m][i]=iagemax+2; |
if(agev[m][i]==1) agev[m][i]=iagemax+2; |
| if((int)agev[m][i] <iagemin-AGEMARGE || (int)agev[m][i] >iagemax+3+AGEMARGE){ |
if((int)agev[m][i] <iagemin-AGEMARGE || (int)agev[m][i] >iagemax+3+AGEMARGE){ |
| printf("Error on individual # %d agev[m][i]=%f <%d-%d or > %d+3+%d m=%d; either change agemin or agemax or fix data\n",i, agev[m][i],iagemin,AGEMARGE, iagemax,AGEMARGE,m); |
printf("Error on individual # %d agev[m][i]=%f <%d-%d or > %d+3+%d m=%d; either change agemin or agemax or fix data\n",i, agev[m][i],iagemin,AGEMARGE, iagemax,AGEMARGE,m); |
| exit(1); |
exit(1); |
| } |
} |
| if (s[m][i]>0 && s[m][i]<=nlstate) { |
if (s[m][i]>0 && s[m][i]<=nlstate) { |
| /*if(i>4620) printf(" i=%d m=%d s[m][i]=%d (int)agev[m][i]=%d weight[i]=%f prop=%f\n",i,m,s[m][i],(int)agev[m][m],weight[i],prop[s[m][i]][(int)agev[m][i]]);*/ |
/*if(i>4620) printf(" i=%d m=%d s[m][i]=%d (int)agev[m][i]=%d weight[i]=%f prop=%f\n",i,m,s[m][i],(int)agev[m][m],weight[i],prop[s[m][i]][(int)agev[m][i]]);*/ |
| prop[s[m][i]][(int)agev[m][i]] += weight[i];/* At age of beginning of transition, where status is known */ |
prop[s[m][i]][(int)agev[m][i]] += weight[i];/* At age of beginning of transition, where status is known */ |
| prop[s[m][i]][iagemax+3] += weight[i]; |
prop[s[m][i]][iagemax+3] += weight[i]; |
| } /* end valid statuses */ |
} /* end valid statuses */ |
| } /* end selection of dates */ |
} /* end selection of dates */ |
| } /* end selection of waves */ |
} /* end selection of waves */ |
| } /* end effective waves */ |
} /* end effective waves */ |
| } /* end bool */ |
} /* end bool */ |
| } |
} |
| for(i=iagemin; i <= iagemax+3; i++){ |
for(i=iagemin; i <= iagemax+3; i++){ |
| for(jk=1,posprop=0; jk <=nlstate ; jk++) { |
for(jk=1,posprop=0; jk <=nlstate ; jk++) { |
| posprop += prop[jk][i]; |
posprop += prop[jk][i]; |
| } |
} |
| |
|
| for(jk=1; jk <=nlstate ; jk++){ |
for(jk=1; jk <=nlstate ; jk++){ |
| if( i <= iagemax){ |
if( i <= iagemax){ |
| if(posprop>=1.e-5){ |
if(posprop>=1.e-5){ |
| probs[i][jk][j1]= prop[jk][i]/posprop; |
probs[i][jk][j1]= prop[jk][i]/posprop; |
| } else{ |
} else{ |
| if(first==1){ |
if(first==1){ |
| first=0; |
first=0; |
| printf("Warning Observed prevalence probs[%d][%d][%d]=%lf because of lack of cases\nSee others on log file...\n",jk,i,j1,probs[i][jk][j1]); |
printf("Warning Observed prevalence probs[%d][%d][%d]=%lf because of lack of cases\nSee others on log file...\n",jk,i,j1,probs[i][jk][j1]); |
| } |
} |
| } |
} |
| } |
} |
| }/* end jk */ |
}/* end jk */ |
| }/* end i */ |
}/* end i */ |
| /*} *//* end i1 */ |
/*} *//* end i1 */ |
| } /* end j1 */ |
} /* end j1 */ |
| |
|
| /* free_ma3x(freq,-1,nlstate+ndeath,-1,nlstate+ndeath, iagemin, iagemax+3);*/ |
/* free_ma3x(freq,-1,nlstate+ndeath,-1,nlstate+ndeath, iagemin, iagemax+3);*/ |
| /*free_vector(pp,1,nlstate);*/ |
/*free_vector(pp,1,nlstate);*/ |
| free_matrix(prop,1,nlstate, iagemin-AGEMARGE,iagemax+3+AGEMARGE); |
free_matrix(prop,1,nlstate, iagemin-AGEMARGE,iagemax+3+AGEMARGE); |
| } /* End of prevalence */ |
} /* End of prevalence */ |
| |
|
| /************* Waves Concatenation ***************/ |
/************* Waves Concatenation ***************/ |
| |
|
|
Line 4520 void cvevsij(double ***eij, double x[],
|
Line 4523 void cvevsij(double ***eij, double x[],
|
| |
|
| { |
{ |
| /* Covariances of health expectancies eij and of total life expectancies according |
/* Covariances of health expectancies eij and of total life expectancies according |
| to initial status i, ei. . |
to initial status i, ei. . |
| */ |
*/ |
| int i, j, nhstepm, hstepm, h, nstepm, k, cptj, cptj2, i2, j2, ij, ji; |
int i, j, nhstepm, hstepm, h, nstepm, k, cptj, cptj2, i2, j2, ij, ji; |
| int nhstepma, nstepma; /* Decreasing with age */ |
int nhstepma, nstepma; /* Decreasing with age */ |
|
Line 4626 void cvevsij(double ***eij, double x[],
|
Line 4629 void cvevsij(double ***eij, double x[],
|
| decrease memory allocation */ |
decrease memory allocation */ |
| for(theta=1; theta <=npar; theta++){ |
for(theta=1; theta <=npar; theta++){ |
| for(i=1; i<=npar; i++){ |
for(i=1; i<=npar; i++){ |
| xp[i] = x[i] + (i==theta ?delti[theta]:0); |
xp[i] = x[i] + (i==theta ?delti[theta]:0); |
| xm[i] = x[i] - (i==theta ?delti[theta]:0); |
xm[i] = x[i] - (i==theta ?delti[theta]:0); |
| } |
} |
| hpxij(p3matp,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, cij); |
hpxij(p3matp,nhstepm,age,hstepm,xp,nlstate,stepm,oldm,savm, cij); |
| hpxij(p3matm,nhstepm,age,hstepm,xm,nlstate,stepm,oldm,savm, cij); |
hpxij(p3matm,nhstepm,age,hstepm,xm,nlstate,stepm,oldm,savm, cij); |
| |
|
| for(j=1; j<= nlstate; j++){ |
for(j=1; j<= nlstate; j++){ |
| for(i=1; i<=nlstate; i++){ |
for(i=1; i<=nlstate; i++){ |
| for(h=0; h<=nhstepm-1; h++){ |
for(h=0; h<=nhstepm-1; h++){ |
| gp[h][(j-1)*nlstate + i] = (p3matp[i][j][h]+p3matp[i][j][h+1])/2.; |
gp[h][(j-1)*nlstate + i] = (p3matp[i][j][h]+p3matp[i][j][h+1])/2.; |
| gm[h][(j-1)*nlstate + i] = (p3matm[i][j][h]+p3matm[i][j][h+1])/2.; |
gm[h][(j-1)*nlstate + i] = (p3matm[i][j][h]+p3matm[i][j][h+1])/2.; |
| } |
} |
| } |
} |
| } |
} |
| |
|
| for(ij=1; ij<= nlstate*nlstate; ij++) |
for(ij=1; ij<= nlstate*nlstate; ij++) |
| for(h=0; h<=nhstepm-1; h++){ |
for(h=0; h<=nhstepm-1; h++){ |
| gradg[h][theta][ij]= (gp[h][ij]-gm[h][ij])/2./delti[theta]; |
gradg[h][theta][ij]= (gp[h][ij]-gm[h][ij])/2./delti[theta]; |
| } |
} |
| }/* End theta */ |
}/* End theta */ |
| |
|
| |
|
| for(h=0; h<=nhstepm-1; h++) |
for(h=0; h<=nhstepm-1; h++) |
| for(j=1; j<=nlstate*nlstate;j++) |
for(j=1; j<=nlstate*nlstate;j++) |
| for(theta=1; theta <=npar; theta++) |
for(theta=1; theta <=npar; theta++) |
| trgradg[h][j][theta]=gradg[h][theta][j]; |
trgradg[h][j][theta]=gradg[h][theta][j]; |
| |
|
| |
|
| for(ij=1;ij<=nlstate*nlstate;ij++) |
for(ij=1;ij<=nlstate*nlstate;ij++) |
| for(ji=1;ji<=nlstate*nlstate;ji++) |
for(ji=1;ji<=nlstate*nlstate;ji++) |
| varhe[ij][ji][(int)age] =0.; |
varhe[ij][ji][(int)age] =0.; |
| |
|
| printf("%d|",(int)age);fflush(stdout); |
printf("%d|",(int)age);fflush(stdout); |
| fprintf(ficlog,"%d|",(int)age);fflush(ficlog); |
fprintf(ficlog,"%d|",(int)age);fflush(ficlog); |
| for(h=0;h<=nhstepm-1;h++){ |
for(h=0;h<=nhstepm-1;h++){ |
| for(k=0;k<=nhstepm-1;k++){ |
for(k=0;k<=nhstepm-1;k++){ |
| matprod2(dnewm,trgradg[h],1,nlstate*nlstate,1,npar,1,npar,matcov); |
matprod2(dnewm,trgradg[h],1,nlstate*nlstate,1,npar,1,npar,matcov); |
| matprod2(doldm,dnewm,1,nlstate*nlstate,1,npar,1,nlstate*nlstate,gradg[k]); |
matprod2(doldm,dnewm,1,nlstate*nlstate,1,npar,1,nlstate*nlstate,gradg[k]); |
| for(ij=1;ij<=nlstate*nlstate;ij++) |
for(ij=1;ij<=nlstate*nlstate;ij++) |
| for(ji=1;ji<=nlstate*nlstate;ji++) |
for(ji=1;ji<=nlstate*nlstate;ji++) |
| varhe[ij][ji][(int)age] += doldm[ij][ji]*hf*hf; |
varhe[ij][ji][(int)age] += doldm[ij][ji]*hf*hf; |
| } |
} |
| } |
} |
| |
|
|
Line 4674 void cvevsij(double ***eij, double x[],
|
Line 4677 void cvevsij(double ***eij, double x[],
|
| hpxij(p3matm,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, cij); |
hpxij(p3matm,nhstepm,age,hstepm,x,nlstate,stepm,oldm, savm, cij); |
| for(i=1; i<=nlstate;i++) |
for(i=1; i<=nlstate;i++) |
| for(j=1; j<=nlstate;j++) |
for(j=1; j<=nlstate;j++) |
| for (h=0, eij[i][j][(int)age]=0; h<=nhstepm-1; h++){ |
for (h=0, eij[i][j][(int)age]=0; h<=nhstepm-1; h++){ |
| eij[i][j][(int)age] += (p3matm[i][j][h]+p3matm[i][j][h+1])/2.0*hf; |
eij[i][j][(int)age] += (p3matm[i][j][h]+p3matm[i][j][h+1])/2.0*hf; |
| |
|
| /* if((int)age==70)printf("i=%2d,j=%2d,h=%2d,age=%3d,%9.4f,%9.4f,%9.4f\n",i,j,h,(int)age,p3mat[i][j][h],hf,eij[i][j][(int)age]);*/ |
/* if((int)age==70)printf("i=%2d,j=%2d,h=%2d,age=%3d,%9.4f,%9.4f,%9.4f\n",i,j,h,(int)age,p3mat[i][j][h],hf,eij[i][j][(int)age]);*/ |
| |
|
| } |
} |
| |
|
| fprintf(ficresstdeij,"%3.0f",age ); |
fprintf(ficresstdeij,"%3.0f",age ); |
| for(i=1; i<=nlstate;i++){ |
for(i=1; i<=nlstate;i++){ |
| eip=0.; |
eip=0.; |
| vip=0.; |
vip=0.; |
| for(j=1; j<=nlstate;j++){ |
for(j=1; j<=nlstate;j++){ |
| eip += eij[i][j][(int)age]; |
eip += eij[i][j][(int)age]; |
| for(k=1; k<=nlstate;k++) /* Sum on j and k of cov(eij,eik) */ |
for(k=1; k<=nlstate;k++) /* Sum on j and k of cov(eij,eik) */ |
| vip += varhe[(j-1)*nlstate+i][(k-1)*nlstate+i][(int)age]; |
vip += varhe[(j-1)*nlstate+i][(k-1)*nlstate+i][(int)age]; |
| fprintf(ficresstdeij," %9.4f (%.4f)", eij[i][j][(int)age], sqrt(varhe[(j-1)*nlstate+i][(j-1)*nlstate+i][(int)age]) ); |
fprintf(ficresstdeij," %9.4f (%.4f)", eij[i][j][(int)age], sqrt(varhe[(j-1)*nlstate+i][(j-1)*nlstate+i][(int)age]) ); |
| } |
} |
| fprintf(ficresstdeij," %9.4f (%.4f)", eip, sqrt(vip)); |
fprintf(ficresstdeij," %9.4f (%.4f)", eip, sqrt(vip)); |
| } |
} |
|
Line 4698 void cvevsij(double ***eij, double x[],
|
Line 4701 void cvevsij(double ***eij, double x[],
|
| fprintf(ficrescveij,"%3.0f",age ); |
fprintf(ficrescveij,"%3.0f",age ); |
| for(i=1; i<=nlstate;i++) |
for(i=1; i<=nlstate;i++) |
| for(j=1; j<=nlstate;j++){ |
for(j=1; j<=nlstate;j++){ |
| cptj= (j-1)*nlstate+i; |
cptj= (j-1)*nlstate+i; |
| for(i2=1; i2<=nlstate;i2++) |
for(i2=1; i2<=nlstate;i2++) |
| for(j2=1; j2<=nlstate;j2++){ |
for(j2=1; j2<=nlstate;j2++){ |
| cptj2= (j2-1)*nlstate+i2; |
cptj2= (j2-1)*nlstate+i2; |
| if(cptj2 <= cptj) |
if(cptj2 <= cptj) |
| fprintf(ficrescveij," %.4f", varhe[cptj][cptj2][(int)age]); |
fprintf(ficrescveij," %.4f", varhe[cptj][cptj2][(int)age]); |
| } |
} |
| } |
} |
| fprintf(ficrescveij,"\n"); |
fprintf(ficrescveij,"\n"); |
| |
|
|
Line 5161 void cvevsij(double ***eij, double x[],
|
Line 5164 void cvevsij(double ***eij, double x[],
|
| |
|
| /************ Variance of one-step probabilities ******************/ |
/************ Variance of one-step probabilities ******************/ |
| void varprob(char optionfilefiname[], double **matcov, double x[], double delti[], int nlstate, double bage, double fage, int ij, int *Tvar, int **nbcode, int *ncodemax, char strstart[]) |
void varprob(char optionfilefiname[], double **matcov, double x[], double delti[], int nlstate, double bage, double fage, int ij, int *Tvar, int **nbcode, int *ncodemax, char strstart[]) |
| { |
{ |
| int i, j=0, k1, l1, tj; |
int i, j=0, k1, l1, tj; |
| int k2, l2, j1, z1; |
int k2, l2, j1, z1; |
| int k=0, l; |
int k=0, l; |
| int first=1, first1, first2; |
int first=1, first1, first2; |
| double cv12, mu1, mu2, lc1, lc2, v12, v21, v11, v22,v1,v2, c12, tnalp; |
double cv12, mu1, mu2, lc1, lc2, v12, v21, v11, v22,v1,v2, c12, tnalp; |
| double **dnewm,**doldm; |
double **dnewm,**doldm; |
| double *xp; |
double *xp; |
| double *gp, *gm; |
double *gp, *gm; |
| double **gradg, **trgradg; |
double **gradg, **trgradg; |
| double **mu; |
double **mu; |
| double age, cov[NCOVMAX+1]; |
double age, cov[NCOVMAX+1]; |
| double std=2.0; /* Number of standard deviation wide of confidence ellipsoids */ |
double std=2.0; /* Number of standard deviation wide of confidence ellipsoids */ |
| int theta; |
int theta; |
| char fileresprob[FILENAMELENGTH]; |
char fileresprob[FILENAMELENGTH]; |
| char fileresprobcov[FILENAMELENGTH]; |
char fileresprobcov[FILENAMELENGTH]; |
| char fileresprobcor[FILENAMELENGTH]; |
char fileresprobcor[FILENAMELENGTH]; |
| double ***varpij; |
double ***varpij; |
| |
|
| strcpy(fileresprob,"PROB_"); |
strcpy(fileresprob,"PROB_"); |
| strcat(fileresprob,fileres); |
strcat(fileresprob,fileres); |
| if((ficresprob=fopen(fileresprob,"w"))==NULL) { |
if((ficresprob=fopen(fileresprob,"w"))==NULL) { |
| printf("Problem with resultfile: %s\n", fileresprob); |
printf("Problem with resultfile: %s\n", fileresprob); |
| fprintf(ficlog,"Problem with resultfile: %s\n", fileresprob); |
fprintf(ficlog,"Problem with resultfile: %s\n", fileresprob); |
| } |
} |
| strcpy(fileresprobcov,"PROBCOV_"); |
strcpy(fileresprobcov,"PROBCOV_"); |
| strcat(fileresprobcov,fileresu); |
strcat(fileresprobcov,fileresu); |
| if((ficresprobcov=fopen(fileresprobcov,"w"))==NULL) { |
if((ficresprobcov=fopen(fileresprobcov,"w"))==NULL) { |
| printf("Problem with resultfile: %s\n", fileresprobcov); |
printf("Problem with resultfile: %s\n", fileresprobcov); |
| fprintf(ficlog,"Problem with resultfile: %s\n", fileresprobcov); |
fprintf(ficlog,"Problem with resultfile: %s\n", fileresprobcov); |
| } |
} |
| strcpy(fileresprobcor,"PROBCOR_"); |
strcpy(fileresprobcor,"PROBCOR_"); |
| strcat(fileresprobcor,fileresu); |
strcat(fileresprobcor,fileresu); |
| if((ficresprobcor=fopen(fileresprobcor,"w"))==NULL) { |
if((ficresprobcor=fopen(fileresprobcor,"w"))==NULL) { |
| printf("Problem with resultfile: %s\n", fileresprobcor); |
printf("Problem with resultfile: %s\n", fileresprobcor); |
| fprintf(ficlog,"Problem with resultfile: %s\n", fileresprobcor); |
fprintf(ficlog,"Problem with resultfile: %s\n", fileresprobcor); |
| } |
} |
| printf("Computing standard deviation of one-step probabilities: result on file '%s' \n",fileresprob); |
printf("Computing standard deviation of one-step probabilities: result on file '%s' \n",fileresprob); |
| fprintf(ficlog,"Computing standard deviation of one-step probabilities: result on file '%s' \n",fileresprob); |
fprintf(ficlog,"Computing standard deviation of one-step probabilities: result on file '%s' \n",fileresprob); |
| printf("Computing matrix of variance covariance of one-step probabilities: result on file '%s' \n",fileresprobcov); |
printf("Computing matrix of variance covariance of one-step probabilities: result on file '%s' \n",fileresprobcov); |
| fprintf(ficlog,"Computing matrix of variance covariance of one-step probabilities: result on file '%s' \n",fileresprobcov); |
fprintf(ficlog,"Computing matrix of variance covariance of one-step probabilities: result on file '%s' \n",fileresprobcov); |
| printf("and correlation matrix of one-step probabilities: result on file '%s' \n",fileresprobcor); |
printf("and correlation matrix of one-step probabilities: result on file '%s' \n",fileresprobcor); |
| fprintf(ficlog,"and correlation matrix of one-step probabilities: result on file '%s' \n",fileresprobcor); |
fprintf(ficlog,"and correlation matrix of one-step probabilities: result on file '%s' \n",fileresprobcor); |
| pstamp(ficresprob); |
pstamp(ficresprob); |
| fprintf(ficresprob,"#One-step probabilities and stand. devi in ()\n"); |
fprintf(ficresprob,"#One-step probabilities and stand. devi in ()\n"); |
| fprintf(ficresprob,"# Age"); |
fprintf(ficresprob,"# Age"); |
| pstamp(ficresprobcov); |
pstamp(ficresprobcov); |
| fprintf(ficresprobcov,"#One-step probabilities and covariance matrix\n"); |
fprintf(ficresprobcov,"#One-step probabilities and covariance matrix\n"); |
| fprintf(ficresprobcov,"# Age"); |
fprintf(ficresprobcov,"# Age"); |
| pstamp(ficresprobcor); |
pstamp(ficresprobcor); |
| fprintf(ficresprobcor,"#One-step probabilities and correlation matrix\n"); |
fprintf(ficresprobcor,"#One-step probabilities and correlation matrix\n"); |
| fprintf(ficresprobcor,"# Age"); |
fprintf(ficresprobcor,"# Age"); |
| |
|
| |
|
| for(i=1; i<=nlstate;i++) |
|
| for(j=1; j<=(nlstate+ndeath);j++){ |
|
| fprintf(ficresprob," p%1d-%1d (SE)",i,j); |
|
| fprintf(ficresprobcov," p%1d-%1d ",i,j); |
|
| fprintf(ficresprobcor," p%1d-%1d ",i,j); |
|
| } |
|
| /* fprintf(ficresprob,"\n"); |
|
| fprintf(ficresprobcov,"\n"); |
|
| fprintf(ficresprobcor,"\n"); |
|
| */ |
|
| xp=vector(1,npar); |
|
| dnewm=matrix(1,(nlstate)*(nlstate+ndeath),1,npar); |
|
| doldm=matrix(1,(nlstate)*(nlstate+ndeath),1,(nlstate)*(nlstate+ndeath)); |
|
| mu=matrix(1,(nlstate)*(nlstate+ndeath), (int) bage, (int)fage); |
|
| varpij=ma3x(1,nlstate*(nlstate+ndeath),1,nlstate*(nlstate+ndeath),(int) bage, (int) fage); |
|
| first=1; |
|
| fprintf(ficgp,"\n# Routine varprob"); |
|
| fprintf(fichtm,"\n<li><h4> Computing and drawing one step probabilities with their confidence intervals</h4></li>\n"); |
|
| fprintf(fichtm,"\n"); |
|
| |
|
| fprintf(fichtm,"\n<li><h4> <a href=\"%s\">Matrix of variance-covariance of one-step probabilities (drawings)</a></h4> this page is important in order to visualize confidence intervals and especially correlation between disability and recovery, or more generally, way in and way back.</li>\n",optionfilehtmcov); |
for(i=1; i<=nlstate;i++) |
| fprintf(fichtmcov,"Current page is file <a href=\"%s\">%s</a><br>\n\n<h4>Matrix of variance-covariance of pairs of step probabilities</h4>\n",optionfilehtmcov, optionfilehtmcov); |
for(j=1; j<=(nlstate+ndeath);j++){ |
| fprintf(fichtmcov,"\nEllipsoids of confidence centered on point (p<inf>ij</inf>, p<inf>kl</inf>) are estimated \ |
fprintf(ficresprob," p%1d-%1d (SE)",i,j); |
| |
fprintf(ficresprobcov," p%1d-%1d ",i,j); |
| |
fprintf(ficresprobcor," p%1d-%1d ",i,j); |
| |
} |
| |
/* fprintf(ficresprob,"\n"); |
| |
fprintf(ficresprobcov,"\n"); |
| |
fprintf(ficresprobcor,"\n"); |
| |
*/ |
| |
xp=vector(1,npar); |
| |
dnewm=matrix(1,(nlstate)*(nlstate+ndeath),1,npar); |
| |
doldm=matrix(1,(nlstate)*(nlstate+ndeath),1,(nlstate)*(nlstate+ndeath)); |
| |
mu=matrix(1,(nlstate)*(nlstate+ndeath), (int) bage, (int)fage); |
| |
varpij=ma3x(1,nlstate*(nlstate+ndeath),1,nlstate*(nlstate+ndeath),(int) bage, (int) fage); |
| |
first=1; |
| |
fprintf(ficgp,"\n# Routine varprob"); |
| |
fprintf(fichtm,"\n<li><h4> Computing and drawing one step probabilities with their confidence intervals</h4></li>\n"); |
| |
fprintf(fichtm,"\n"); |
| |
|
| |
fprintf(fichtm,"\n<li><h4> <a href=\"%s\">Matrix of variance-covariance of one-step probabilities (drawings)</a></h4> this page is important in order to visualize confidence intervals and especially correlation between disability and recovery, or more generally, way in and way back.</li>\n",optionfilehtmcov); |
| |
fprintf(fichtmcov,"Current page is file <a href=\"%s\">%s</a><br>\n\n<h4>Matrix of variance-covariance of pairs of step probabilities</h4>\n",optionfilehtmcov, optionfilehtmcov); |
| |
fprintf(fichtmcov,"\nEllipsoids of confidence centered on point (p<inf>ij</inf>, p<inf>kl</inf>) are estimated \ |
| and drawn. It helps understanding how is the covariance between two incidences.\ |
and drawn. It helps understanding how is the covariance between two incidences.\ |
| They are expressed in year<sup>-1</sup> in order to be less dependent of stepm.<br>\n"); |
They are expressed in year<sup>-1</sup> in order to be less dependent of stepm.<br>\n"); |
| fprintf(fichtmcov,"\n<br> Contour plot corresponding to x'cov<sup>-1</sup>x = 4 (where x is the column vector (pij,pkl)) are drawn. \ |
fprintf(fichtmcov,"\n<br> Contour plot corresponding to x'cov<sup>-1</sup>x = 4 (where x is the column vector (pij,pkl)) are drawn. \ |
| It can be understood this way: if pij and pkl where uncorrelated the (2x2) matrix of covariance \ |
It can be understood this way: if pij and pkl where uncorrelated the (2x2) matrix of covariance \ |
| would have been (1/(var pij), 0 , 0, 1/(var pkl)), and the confidence interval would be 2 \ |
would have been (1/(var pij), 0 , 0, 1/(var pkl)), and the confidence interval would be 2 \ |
| standard deviations wide on each axis. <br>\ |
standard deviations wide on each axis. <br>\ |
|
Line 5248 standard deviations wide on each axis. <
|
Line 5251 standard deviations wide on each axis. <
|
| and made the appropriate rotation to look at the uncorrelated principal directions.<br>\ |
and made the appropriate rotation to look at the uncorrelated principal directions.<br>\ |
| To be simple, these graphs help to understand the significativity of each parameter in relation to a second other one.<br> \n"); |
To be simple, these graphs help to understand the significativity of each parameter in relation to a second other one.<br> \n"); |
| |
|
| cov[1]=1; |
cov[1]=1; |
| /* tj=cptcoveff; */ |
/* tj=cptcoveff; */ |
| tj = (int) pow(2,cptcoveff); |
tj = (int) pow(2,cptcoveff); |
| if (cptcovn<1) {tj=1;ncodemax[1]=1;} |
if (cptcovn<1) {tj=1;ncodemax[1]=1;} |
| j1=0; |
j1=0; |
| for(j1=1; j1<=tj;j1++){ /* For each valid combination of covariates */ |
for(j1=1; j1<=tj;j1++){ /* For each valid combination of covariates */ |
| if (cptcovn>0) { |
if (cptcovn>0) { |
| fprintf(ficresprob, "\n#********** Variable "); |
fprintf(ficresprob, "\n#********** Variable "); |
| for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprob, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprob, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
| fprintf(ficresprob, "**********\n#\n"); |
fprintf(ficresprob, "**********\n#\n"); |
| fprintf(ficresprobcov, "\n#********** Variable "); |
fprintf(ficresprobcov, "\n#********** Variable "); |
| for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprobcov, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprobcov, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
| fprintf(ficresprobcov, "**********\n#\n"); |
fprintf(ficresprobcov, "**********\n#\n"); |
| |
|
| fprintf(ficgp, "\n#********** Variable "); |
fprintf(ficgp, "\n#********** Variable "); |
| for (z1=1; z1<=cptcoveff; z1++) fprintf(ficgp, " V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
for (z1=1; z1<=cptcoveff; z1++) fprintf(ficgp, " V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
| fprintf(ficgp, "**********\n#\n"); |
fprintf(ficgp, "**********\n#\n"); |
| |
|
| |
|
| fprintf(fichtmcov, "\n<hr size=\"2\" color=\"#EC5E5E\">********** Variable "); |
fprintf(fichtmcov, "\n<hr size=\"2\" color=\"#EC5E5E\">********** Variable "); |
| for (z1=1; z1<=cptcoveff; z1++) fprintf(fichtm, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
for (z1=1; z1<=cptcoveff; z1++) fprintf(fichtm, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
| fprintf(fichtmcov, "**********\n<hr size=\"2\" color=\"#EC5E5E\">"); |
fprintf(fichtmcov, "**********\n<hr size=\"2\" color=\"#EC5E5E\">"); |
| |
|
| fprintf(ficresprobcor, "\n#********** Variable "); |
fprintf(ficresprobcor, "\n#********** Variable "); |
| for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprobcor, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresprobcor, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtabm(j1,z1)]); |
| fprintf(ficresprobcor, "**********\n#"); |
fprintf(ficresprobcor, "**********\n#"); |
| if(invalidvarcomb[j1]){ |
if(invalidvarcomb[j1]){ |
| fprintf(ficgp,"\n#Combination (%d) ignored because no cases \n",j1); |
fprintf(ficgp,"\n#Combination (%d) ignored because no cases \n",j1); |
| fprintf(fichtmcov,"\n<h3>Combination (%d) ignored because no cases </h3>\n",j1); |
fprintf(fichtmcov,"\n<h3>Combination (%d) ignored because no cases </h3>\n",j1); |
| continue; |
continue; |
| } |
} |
| } |
} |
| gradg=matrix(1,npar,1,(nlstate)*(nlstate+ndeath)); |
gradg=matrix(1,npar,1,(nlstate)*(nlstate+ndeath)); |
| trgradg=matrix(1,(nlstate)*(nlstate+ndeath),1,npar); |
trgradg=matrix(1,(nlstate)*(nlstate+ndeath),1,npar); |
| gp=vector(1,(nlstate)*(nlstate+ndeath)); |
gp=vector(1,(nlstate)*(nlstate+ndeath)); |
| gm=vector(1,(nlstate)*(nlstate+ndeath)); |
gm=vector(1,(nlstate)*(nlstate+ndeath)); |
| for (age=bage; age<=fage; age ++){ |
for (age=bage; age<=fage; age ++){ |
| cov[2]=age; |
cov[2]=age; |
| if(nagesqr==1) |
if(nagesqr==1) |
| cov[3]= age*age; |
cov[3]= age*age; |
| for (k=1; k<=cptcovn;k++) { |
for (k=1; k<=cptcovn;k++) { |
| cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,k)]; |
cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,k)]; |
| /*cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,Tvar[k])];*//* j1 1 2 3 4 |
/*cov[2+nagesqr+k]=nbcode[Tvar[k]][codtabm(j1,Tvar[k])];*//* j1 1 2 3 4 |
| * 1 1 1 1 1 |
* 1 1 1 1 1 |
| * 2 2 1 1 1 |
* 2 2 1 1 1 |
| * 3 1 2 1 1 |
* 3 1 2 1 1 |
| */ |
*/ |
| /* nbcode[1][1]=0 nbcode[1][2]=1;*/ |
/* nbcode[1][1]=0 nbcode[1][2]=1;*/ |
| } |
} |
| /* for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */ |
/* for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2]; */ |
| for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; |
for (k=1; k<=cptcovage;k++) cov[2+Tage[k]]=nbcode[Tvar[Tage[k]]][codtabm(ij,k)]*cov[2]; |
| for (k=1; k<=cptcovprod;k++) |
for (k=1; k<=cptcovprod;k++) |
| cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)]; |
cov[2+nagesqr+Tprod[k]]=nbcode[Tvard[k][1]][codtabm(ij,k)]*nbcode[Tvard[k][2]][codtabm(ij,k)]; |
| |
|
| |
|
| for(theta=1; theta <=npar; theta++){ |
for(theta=1; theta <=npar; theta++){ |
| for(i=1; i<=npar; i++) |
for(i=1; i<=npar; i++) |
| xp[i] = x[i] + (i==theta ?delti[theta]:(double)0); |
xp[i] = x[i] + (i==theta ?delti[theta]:(double)0); |
| |
|
| pmij(pmmij,cov,ncovmodel,xp,nlstate); |
pmij(pmmij,cov,ncovmodel,xp,nlstate); |
| |
|
| k=0; |
k=0; |
| for(i=1; i<= (nlstate); i++){ |
for(i=1; i<= (nlstate); i++){ |
| for(j=1; j<=(nlstate+ndeath);j++){ |
for(j=1; j<=(nlstate+ndeath);j++){ |
| k=k+1; |
k=k+1; |
| gp[k]=pmmij[i][j]; |
gp[k]=pmmij[i][j]; |
| } |
} |
| } |
} |
| |
|
| for(i=1; i<=npar; i++) |
for(i=1; i<=npar; i++) |
| xp[i] = x[i] - (i==theta ?delti[theta]:(double)0); |
xp[i] = x[i] - (i==theta ?delti[theta]:(double)0); |
| |
|
| pmij(pmmij,cov,ncovmodel,xp,nlstate); |
pmij(pmmij,cov,ncovmodel,xp,nlstate); |
| k=0; |
k=0; |
| for(i=1; i<=(nlstate); i++){ |
for(i=1; i<=(nlstate); i++){ |
| for(j=1; j<=(nlstate+ndeath);j++){ |
for(j=1; j<=(nlstate+ndeath);j++){ |
| k=k+1; |
k=k+1; |
| gm[k]=pmmij[i][j]; |
gm[k]=pmmij[i][j]; |
| } |
} |
| } |
} |
| |
|
| for(i=1; i<= (nlstate)*(nlstate+ndeath); i++) |
for(i=1; i<= (nlstate)*(nlstate+ndeath); i++) |
| gradg[theta][i]=(gp[i]-gm[i])/(double)2./delti[theta]; |
gradg[theta][i]=(gp[i]-gm[i])/(double)2./delti[theta]; |
| } |
} |
| |
|
| for(j=1; j<=(nlstate)*(nlstate+ndeath);j++) |
for(j=1; j<=(nlstate)*(nlstate+ndeath);j++) |
| for(theta=1; theta <=npar; theta++) |
for(theta=1; theta <=npar; theta++) |
| trgradg[j][theta]=gradg[theta][j]; |
trgradg[j][theta]=gradg[theta][j]; |
| |
|
| matprod2(dnewm,trgradg,1,(nlstate)*(nlstate+ndeath),1,npar,1,npar,matcov); |
matprod2(dnewm,trgradg,1,(nlstate)*(nlstate+ndeath),1,npar,1,npar,matcov); |
| matprod2(doldm,dnewm,1,(nlstate)*(nlstate+ndeath),1,npar,1,(nlstate)*(nlstate+ndeath),gradg); |
matprod2(doldm,dnewm,1,(nlstate)*(nlstate+ndeath),1,npar,1,(nlstate)*(nlstate+ndeath),gradg); |
| |
|
| pmij(pmmij,cov,ncovmodel,x,nlstate); |
pmij(pmmij,cov,ncovmodel,x,nlstate); |
| |
|
| k=0; |
k=0; |
| for(i=1; i<=(nlstate); i++){ |
for(i=1; i<=(nlstate); i++){ |
| for(j=1; j<=(nlstate+ndeath);j++){ |
for(j=1; j<=(nlstate+ndeath);j++){ |
| k=k+1; |
k=k+1; |
| mu[k][(int) age]=pmmij[i][j]; |
mu[k][(int) age]=pmmij[i][j]; |
| } |
} |
| } |
} |
| for(i=1;i<=(nlstate)*(nlstate+ndeath);i++) |
for(i=1;i<=(nlstate)*(nlstate+ndeath);i++) |
| for(j=1;j<=(nlstate)*(nlstate+ndeath);j++) |
for(j=1;j<=(nlstate)*(nlstate+ndeath);j++) |
| varpij[i][j][(int)age] = doldm[i][j]; |
varpij[i][j][(int)age] = doldm[i][j]; |
| |
|
| /*printf("\n%d ",(int)age); |
/*printf("\n%d ",(int)age); |
| for (i=1; i<=(nlstate)*(nlstate+ndeath);i++){ |
for (i=1; i<=(nlstate)*(nlstate+ndeath);i++){ |
| printf("%e [%e ;%e] ",gm[i],gm[i]-2*sqrt(doldm[i][i]),gm[i]+2*sqrt(doldm[i][i])); |
printf("%e [%e ;%e] ",gm[i],gm[i]-2*sqrt(doldm[i][i]),gm[i]+2*sqrt(doldm[i][i])); |
| fprintf(ficlog,"%e [%e ;%e] ",gm[i],gm[i]-2*sqrt(doldm[i][i]),gm[i]+2*sqrt(doldm[i][i])); |
fprintf(ficlog,"%e [%e ;%e] ",gm[i],gm[i]-2*sqrt(doldm[i][i]),gm[i]+2*sqrt(doldm[i][i])); |
| }*/ |
}*/ |
| |
|
| fprintf(ficresprob,"\n%d ",(int)age); |
fprintf(ficresprob,"\n%d ",(int)age); |
| fprintf(ficresprobcov,"\n%d ",(int)age); |
fprintf(ficresprobcov,"\n%d ",(int)age); |
| fprintf(ficresprobcor,"\n%d ",(int)age); |
fprintf(ficresprobcor,"\n%d ",(int)age); |
| |
|
| for (i=1; i<=(nlstate)*(nlstate+ndeath);i++) |
for (i=1; i<=(nlstate)*(nlstate+ndeath);i++) |
| fprintf(ficresprob,"%11.3e (%11.3e) ",mu[i][(int) age],sqrt(varpij[i][i][(int)age])); |
fprintf(ficresprob,"%11.3e (%11.3e) ",mu[i][(int) age],sqrt(varpij[i][i][(int)age])); |
| for (i=1; i<=(nlstate)*(nlstate+ndeath);i++){ |
for (i=1; i<=(nlstate)*(nlstate+ndeath);i++){ |
| fprintf(ficresprobcov,"%11.3e ",mu[i][(int) age]); |
fprintf(ficresprobcov,"%11.3e ",mu[i][(int) age]); |
| fprintf(ficresprobcor,"%11.3e ",mu[i][(int) age]); |
fprintf(ficresprobcor,"%11.3e ",mu[i][(int) age]); |
| } |
} |
| i=0; |
i=0; |
| for (k=1; k<=(nlstate);k++){ |
for (k=1; k<=(nlstate);k++){ |
| for (l=1; l<=(nlstate+ndeath);l++){ |
for (l=1; l<=(nlstate+ndeath);l++){ |
| i++; |
i++; |
| fprintf(ficresprobcov,"\n%d %d-%d",(int)age,k,l); |
fprintf(ficresprobcov,"\n%d %d-%d",(int)age,k,l); |
| fprintf(ficresprobcor,"\n%d %d-%d",(int)age,k,l); |
fprintf(ficresprobcor,"\n%d %d-%d",(int)age,k,l); |
| for (j=1; j<=i;j++){ |
for (j=1; j<=i;j++){ |
| /* printf(" k=%d l=%d i=%d j=%d\n",k,l,i,j);fflush(stdout); */ |
/* printf(" k=%d l=%d i=%d j=%d\n",k,l,i,j);fflush(stdout); */ |
| fprintf(ficresprobcov," %11.3e",varpij[i][j][(int)age]); |
fprintf(ficresprobcov," %11.3e",varpij[i][j][(int)age]); |
| fprintf(ficresprobcor," %11.3e",varpij[i][j][(int) age]/sqrt(varpij[i][i][(int) age])/sqrt(varpij[j][j][(int)age])); |
fprintf(ficresprobcor," %11.3e",varpij[i][j][(int) age]/sqrt(varpij[i][i][(int) age])/sqrt(varpij[j][j][(int)age])); |
| } |
} |
| } |
} |
| }/* end of loop for state */ |
}/* end of loop for state */ |
| } /* end of loop for age */ |
} /* end of loop for age */ |
| free_vector(gp,1,(nlstate+ndeath)*(nlstate+ndeath)); |
free_vector(gp,1,(nlstate+ndeath)*(nlstate+ndeath)); |
| free_vector(gm,1,(nlstate+ndeath)*(nlstate+ndeath)); |
free_vector(gm,1,(nlstate+ndeath)*(nlstate+ndeath)); |
| free_matrix(trgradg,1,(nlstate+ndeath)*(nlstate+ndeath),1,npar); |
free_matrix(trgradg,1,(nlstate+ndeath)*(nlstate+ndeath),1,npar); |
| free_matrix(gradg,1,(nlstate+ndeath)*(nlstate+ndeath),1,npar); |
free_matrix(gradg,1,(nlstate+ndeath)*(nlstate+ndeath),1,npar); |
| |
|
| /* Confidence intervalle of pij */ |
/* Confidence intervalle of pij */ |
| /* |
/* |
| fprintf(ficgp,"\nunset parametric;unset label"); |
fprintf(ficgp,"\nunset parametric;unset label"); |
| fprintf(ficgp,"\nset log y;unset log x; set xlabel \"Age\";set ylabel \"probability (year-1)\""); |
fprintf(ficgp,"\nset log y;unset log x; set xlabel \"Age\";set ylabel \"probability (year-1)\""); |
| fprintf(ficgp,"\nset ter png small\nset size 0.65,0.65"); |
fprintf(ficgp,"\nset ter png small\nset size 0.65,0.65"); |
| fprintf(fichtm,"\n<br>Probability with confidence intervals expressed in year<sup>-1</sup> :<a href=\"pijgr%s.png\">pijgr%s.png</A>, ",optionfilefiname,optionfilefiname); |
fprintf(fichtm,"\n<br>Probability with confidence intervals expressed in year<sup>-1</sup> :<a href=\"pijgr%s.png\">pijgr%s.png</A>, ",optionfilefiname,optionfilefiname); |
| fprintf(fichtm,"\n<br><img src=\"pijgr%s.png\"> ",optionfilefiname); |
fprintf(fichtm,"\n<br><img src=\"pijgr%s.png\"> ",optionfilefiname); |
| fprintf(ficgp,"\nset out \"pijgr%s.png\"",optionfilefiname); |
fprintf(ficgp,"\nset out \"pijgr%s.png\"",optionfilefiname); |
| fprintf(ficgp,"\nplot \"%s\" every :::%d::%d u 1:2 \"\%%lf",k1,k2,xfilevarprob); |
fprintf(ficgp,"\nplot \"%s\" every :::%d::%d u 1:2 \"\%%lf",k1,k2,xfilevarprob); |
| */ |
*/ |
| |
|
| /* Drawing ellipsoids of confidence of two variables p(k1-l1,k2-l2)*/ |
/* Drawing ellipsoids of confidence of two variables p(k1-l1,k2-l2)*/ |
| first1=1;first2=2; |
first1=1;first2=2; |
| for (k2=1; k2<=(nlstate);k2++){ |
for (k2=1; k2<=(nlstate);k2++){ |
| for (l2=1; l2<=(nlstate+ndeath);l2++){ |
for (l2=1; l2<=(nlstate+ndeath);l2++){ |
| if(l2==k2) continue; |
if(l2==k2) continue; |
| j=(k2-1)*(nlstate+ndeath)+l2; |
j=(k2-1)*(nlstate+ndeath)+l2; |
| for (k1=1; k1<=(nlstate);k1++){ |
for (k1=1; k1<=(nlstate);k1++){ |
| for (l1=1; l1<=(nlstate+ndeath);l1++){ |
for (l1=1; l1<=(nlstate+ndeath);l1++){ |
| if(l1==k1) continue; |
if(l1==k1) continue; |
| i=(k1-1)*(nlstate+ndeath)+l1; |
i=(k1-1)*(nlstate+ndeath)+l1; |
| if(i<=j) continue; |
if(i<=j) continue; |
| for (age=bage; age<=fage; age ++){ |
for (age=bage; age<=fage; age ++){ |
| if ((int)age %5==0){ |
if ((int)age %5==0){ |
| v1=varpij[i][i][(int)age]/stepm*YEARM/stepm*YEARM; |
v1=varpij[i][i][(int)age]/stepm*YEARM/stepm*YEARM; |
| v2=varpij[j][j][(int)age]/stepm*YEARM/stepm*YEARM; |
v2=varpij[j][j][(int)age]/stepm*YEARM/stepm*YEARM; |
| cv12=varpij[i][j][(int)age]/stepm*YEARM/stepm*YEARM; |
cv12=varpij[i][j][(int)age]/stepm*YEARM/stepm*YEARM; |
| mu1=mu[i][(int) age]/stepm*YEARM ; |
mu1=mu[i][(int) age]/stepm*YEARM ; |
| mu2=mu[j][(int) age]/stepm*YEARM; |
mu2=mu[j][(int) age]/stepm*YEARM; |
| c12=cv12/sqrt(v1*v2); |
c12=cv12/sqrt(v1*v2); |
| /* Computing eigen value of matrix of covariance */ |
/* Computing eigen value of matrix of covariance */ |
| lc1=((v1+v2)+sqrt((v1+v2)*(v1+v2) - 4*(v1*v2-cv12*cv12)))/2.; |
lc1=((v1+v2)+sqrt((v1+v2)*(v1+v2) - 4*(v1*v2-cv12*cv12)))/2.; |
| lc2=((v1+v2)-sqrt((v1+v2)*(v1+v2) - 4*(v1*v2-cv12*cv12)))/2.; |
lc2=((v1+v2)-sqrt((v1+v2)*(v1+v2) - 4*(v1*v2-cv12*cv12)))/2.; |
| if ((lc2 <0) || (lc1 <0) ){ |
if ((lc2 <0) || (lc1 <0) ){ |
| if(first2==1){ |
if(first2==1){ |
| first1=0; |
first1=0; |
| printf("Strange: j1=%d One eigen value of 2x2 matrix of covariance is negative, lc1=%11.3e, lc2=%11.3e, v1=%11.3e, v2=%11.3e, cv12=%11.3e.\n It means that the matrix was not well estimated (varpij), for i=%2d, j=%2d, age=%4d .\n See files %s and %s. Probably WRONG RESULTS. See log file for details...\n", j1, lc1, lc2, v1, v2, cv12, i, j, (int)age,fileresprobcov, fileresprobcor); |
printf("Strange: j1=%d One eigen value of 2x2 matrix of covariance is negative, lc1=%11.3e, lc2=%11.3e, v1=%11.3e, v2=%11.3e, cv12=%11.3e.\n It means that the matrix was not well estimated (varpij), for i=%2d, j=%2d, age=%4d .\n See files %s and %s. Probably WRONG RESULTS. See log file for details...\n", j1, lc1, lc2, v1, v2, cv12, i, j, (int)age,fileresprobcov, fileresprobcor); |
| } |
} |
| fprintf(ficlog,"Strange: j1=%d One eigen value of 2x2 matrix of covariance is negative, lc1=%11.3e, lc2=%11.3e, v1=%11.3e, v2=%11.3e, cv12=%11.3e.\n It means that the matrix was not well estimated (varpij), for i=%2d, j=%2d, age=%4d .\n See files %s and %s. Probably WRONG RESULTS.\n", j1, lc1, lc2, v1, v2, cv12, i, j, (int)age,fileresprobcov, fileresprobcor);fflush(ficlog); |
fprintf(ficlog,"Strange: j1=%d One eigen value of 2x2 matrix of covariance is negative, lc1=%11.3e, lc2=%11.3e, v1=%11.3e, v2=%11.3e, cv12=%11.3e.\n It means that the matrix was not well estimated (varpij), for i=%2d, j=%2d, age=%4d .\n See files %s and %s. Probably WRONG RESULTS.\n", j1, lc1, lc2, v1, v2, cv12, i, j, (int)age,fileresprobcov, fileresprobcor);fflush(ficlog); |
| /* lc1=fabs(lc1); */ /* If we want to have them positive */ |
/* lc1=fabs(lc1); */ /* If we want to have them positive */ |
| /* lc2=fabs(lc2); */ |
/* lc2=fabs(lc2); */ |
| } |
} |
| |
|
| /* Eigen vectors */ |
/* Eigen vectors */ |
| v11=(1./sqrt(1+(v1-lc1)*(v1-lc1)/cv12/cv12)); |
v11=(1./sqrt(1+(v1-lc1)*(v1-lc1)/cv12/cv12)); |
| /*v21=sqrt(1.-v11*v11); *//* error */ |
/*v21=sqrt(1.-v11*v11); *//* error */ |
| v21=(lc1-v1)/cv12*v11; |
v21=(lc1-v1)/cv12*v11; |
| v12=-v21; |
v12=-v21; |
| v22=v11; |
v22=v11; |
| tnalp=v21/v11; |
tnalp=v21/v11; |
| if(first1==1){ |
if(first1==1){ |
| first1=0; |
first1=0; |
| printf("%d %d%d-%d%d mu %.4e %.4e Var %.4e %.4e cor %.3f cov %.4e Eig %.3e %.3e 1stv %.3f %.3f tang %.3f\nOthers in log...\n",(int) age,k1,l1,k2,l2,mu1,mu2,v1,v2,c12,cv12,lc1,lc2,v11,v21,tnalp); |
printf("%d %d%d-%d%d mu %.4e %.4e Var %.4e %.4e cor %.3f cov %.4e Eig %.3e %.3e 1stv %.3f %.3f tang %.3f\nOthers in log...\n",(int) age,k1,l1,k2,l2,mu1,mu2,v1,v2,c12,cv12,lc1,lc2,v11,v21,tnalp); |
| } |
} |
| fprintf(ficlog,"%d %d%d-%d%d mu %.4e %.4e Var %.4e %.4e cor %.3f cov %.4e Eig %.3e %.3e 1stv %.3f %.3f tan %.3f\n",(int) age,k1,l1,k2,l2,mu1,mu2,v1,v2,c12,cv12,lc1,lc2,v11,v21,tnalp); |
fprintf(ficlog,"%d %d%d-%d%d mu %.4e %.4e Var %.4e %.4e cor %.3f cov %.4e Eig %.3e %.3e 1stv %.3f %.3f tan %.3f\n",(int) age,k1,l1,k2,l2,mu1,mu2,v1,v2,c12,cv12,lc1,lc2,v11,v21,tnalp); |
| /*printf(fignu*/ |
/*printf(fignu*/ |
| /* mu1+ v11*lc1*cost + v12*lc2*sin(t) */ |
/* mu1+ v11*lc1*cost + v12*lc2*sin(t) */ |
| /* mu2+ v21*lc1*cost + v22*lc2*sin(t) */ |
/* mu2+ v21*lc1*cost + v22*lc2*sin(t) */ |
| if(first==1){ |
if(first==1){ |
| first=0; |
first=0; |
| fprintf(ficgp,"\n# Ellipsoids of confidence\n#\n"); |
fprintf(ficgp,"\n# Ellipsoids of confidence\n#\n"); |
| fprintf(ficgp,"\nset parametric;unset label"); |
fprintf(ficgp,"\nset parametric;unset label"); |
| fprintf(ficgp,"\nset log y;set log x; set xlabel \"p%1d%1d (year-1)\";set ylabel \"p%1d%1d (year-1)\"",k1,l1,k2,l2); |
fprintf(ficgp,"\nset log y;set log x; set xlabel \"p%1d%1d (year-1)\";set ylabel \"p%1d%1d (year-1)\"",k1,l1,k2,l2); |
| fprintf(ficgp,"\nset ter svg size 640, 480"); |
fprintf(ficgp,"\nset ter svg size 640, 480"); |
| fprintf(fichtmcov,"\n<br>Ellipsoids of confidence cov(p%1d%1d,p%1d%1d) expressed in year<sup>-1</sup>\ |
fprintf(fichtmcov,"\n<br>Ellipsoids of confidence cov(p%1d%1d,p%1d%1d) expressed in year<sup>-1</sup>\ |
| :<a href=\"%s_%d%1d%1d-%1d%1d.svg\"> \ |
:<a href=\"%s_%d%1d%1d-%1d%1d.svg\"> \ |
| %s_%d%1d%1d-%1d%1d.svg</A>, ",k1,l1,k2,l2,\ |
%s_%d%1d%1d-%1d%1d.svg</A>, ",k1,l1,k2,l2,\ |
| subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2, \ |
subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2, \ |
| subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
| fprintf(fichtmcov,"\n<br><img src=\"%s_%d%1d%1d-%1d%1d.svg\"> ",subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
fprintf(fichtmcov,"\n<br><img src=\"%s_%d%1d%1d-%1d%1d.svg\"> ",subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
| fprintf(fichtmcov,"\n<br> Correlation at age %d (%.3f),",(int) age, c12); |
fprintf(fichtmcov,"\n<br> Correlation at age %d (%.3f),",(int) age, c12); |
| fprintf(ficgp,"\nset out \"%s_%d%1d%1d-%1d%1d.svg\"",subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
fprintf(ficgp,"\nset out \"%s_%d%1d%1d-%1d%1d.svg\"",subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
| fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu1,mu2); |
fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu1,mu2); |
| fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k1,l1,k2,l2); |
fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k1,l1,k2,l2); |
| fprintf(ficgp,"\nplot [-pi:pi] %11.3e+ %.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)), %11.3e +%.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)) not", \ |
fprintf(ficgp,"\nplot [-pi:pi] %11.3e+ %.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)), %11.3e +%.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)) not", \ |
| mu1,std,v11,sqrt(lc1),v12,sqrt(lc2), \ |
mu1,std,v11,sqrt(lc1),v12,sqrt(lc2), \ |
| mu2,std,v21,sqrt(lc1),v22,sqrt(lc2)); |
mu2,std,v21,sqrt(lc1),v22,sqrt(lc2)); |
| }else{ |
}else{ |
| first=0; |
first=0; |
| fprintf(fichtmcov," %d (%.3f),",(int) age, c12); |
fprintf(fichtmcov," %d (%.3f),",(int) age, c12); |
| fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k1,l1,k2,l2); |
fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k1,l1,k2,l2); |
| fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu1,mu2); |
fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu1,mu2); |
| fprintf(ficgp,"\nreplot %11.3e+ %.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)), %11.3e +%.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)) not", \ |
fprintf(ficgp,"\nreplot %11.3e+ %.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)), %11.3e +%.3f*(%11.3e*%11.3e*cos(t)+%11.3e*%11.3e*sin(t)) not", \ |
| mu1,std,v11,sqrt(lc1),v12,sqrt(lc2), \ |
mu1,std,v11,sqrt(lc1),v12,sqrt(lc2), \ |
| mu2,std,v21,sqrt(lc1),v22,sqrt(lc2)); |
mu2,std,v21,sqrt(lc1),v22,sqrt(lc2)); |
| }/* if first */ |
}/* if first */ |
| } /* age mod 5 */ |
} /* age mod 5 */ |
| } /* end loop age */ |
} /* end loop age */ |
| fprintf(ficgp,"\nset out;\nset out \"%s_%d%1d%1d-%1d%1d.svg\";replot;set out;",subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
fprintf(ficgp,"\nset out;\nset out \"%s_%d%1d%1d-%1d%1d.svg\";replot;set out;",subdirf2(optionfilefiname,"VARPIJGR_"), j1,k1,l1,k2,l2); |
| first=1; |
first=1; |
| } /*l12 */ |
} /*l12 */ |
| } /* k12 */ |
} /* k12 */ |
| } /*l1 */ |
} /*l1 */ |
| }/* k1 */ |
}/* k1 */ |
| } /* loop on combination of covariates j1 */ |
} /* loop on combination of covariates j1 */ |
| free_ma3x(varpij,1,nlstate,1,nlstate+ndeath,(int) bage, (int)fage); |
free_ma3x(varpij,1,nlstate,1,nlstate+ndeath,(int) bage, (int)fage); |
| free_matrix(mu,1,(nlstate+ndeath)*(nlstate+ndeath),(int) bage, (int)fage); |
free_matrix(mu,1,(nlstate+ndeath)*(nlstate+ndeath),(int) bage, (int)fage); |
| free_matrix(doldm,1,(nlstate)*(nlstate+ndeath),1,(nlstate)*(nlstate+ndeath)); |
free_matrix(doldm,1,(nlstate)*(nlstate+ndeath),1,(nlstate)*(nlstate+ndeath)); |
| free_matrix(dnewm,1,(nlstate)*(nlstate+ndeath),1,npar); |
free_matrix(dnewm,1,(nlstate)*(nlstate+ndeath),1,npar); |
| free_vector(xp,1,npar); |
free_vector(xp,1,npar); |
| fclose(ficresprob); |
fclose(ficresprob); |
| fclose(ficresprobcov); |
fclose(ficresprobcov); |
| fclose(ficresprobcor); |
fclose(ficresprobcor); |
| fflush(ficgp); |
fflush(ficgp); |
| fflush(fichtmcov); |
fflush(fichtmcov); |
| } |
} |
| |
|
| |
|
| /******************* Printing html file ***********/ |
/******************* Printing html file ***********/ |
|
Line 5536 void printinghtml(char fileresu[], char
|
Line 5539 void printinghtml(char fileresu[], char
|
| <a href=\"%s\">%s</a> <br>\n</li>", subdirf2(fileresu,"F_"),subdirf2(fileresu,"F_")); |
<a href=\"%s\">%s</a> <br>\n</li>", subdirf2(fileresu,"F_"),subdirf2(fileresu,"F_")); |
| } |
} |
| |
|
| fprintf(fichtm," \n<ul><li><b>Graphs</b></li><p>"); |
fprintf(fichtm," \n<ul><li><b>Graphs</b></li><p>"); |
| |
|
| m=pow(2,cptcoveff); |
m=pow(2,cptcoveff); |
| if (cptcovn < 1) {m=1;ncodemax[1]=1;} |
if (cptcovn < 1) {m=1;ncodemax[1]=1;} |
| |
|
| jj1=0; |
jj1=0; |
| for(k1=1; k1<=m;k1++){ |
for(k1=1; k1<=m;k1++){ |
| |
|
| /* for(i1=1; i1<=ncodemax[k1];i1++){ */ |
/* for(i1=1; i1<=ncodemax[k1];i1++){ */ |
| jj1++; |
jj1++; |
| if (cptcovn > 0) { |
if (cptcovn > 0) { |
| fprintf(fichtm,"<hr size=\"2\" color=\"#EC5E5E\">************ Results for covariates"); |
fprintf(fichtm,"<hr size=\"2\" color=\"#EC5E5E\">************ Results for covariates"); |
| for (cpt=1; cpt<=cptcoveff;cpt++){ |
for (cpt=1; cpt<=cptcoveff;cpt++){ |
| fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]); |
fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]); |
| printf(" V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]);fflush(stdout); |
printf(" V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]);fflush(stdout); |
| } |
} |
| fprintf(fichtm," ************\n<hr size=\"2\" color=\"#EC5E5E\">"); |
fprintf(fichtm," ************\n<hr size=\"2\" color=\"#EC5E5E\">"); |
| if(invalidvarcomb[k1]){ |
if(invalidvarcomb[k1]){ |
| fprintf(fichtm,"\n<h3>Combination (%d) ignored because no cases </h3>\n",k1); |
fprintf(fichtm,"\n<h3>Combination (%d) ignored because no cases </h3>\n",k1); |
| printf("\nCombination (%d) ignored because no cases \n",k1); |
printf("\nCombination (%d) ignored because no cases \n",k1); |
| continue; |
continue; |
| } |
} |
| } |
} |
| /* aij, bij */ |
/* aij, bij */ |
| fprintf(fichtm,"<br>- Logit model (yours is: 1+age+%s), for example: logit(pij)=log(pij/pii)= aij+ bij age + V1 age + etc. as a function of age: <a href=\"%s_%d-1.svg\">%s_%d-1.svg</a><br> \ |
fprintf(fichtm,"<br>- Logit model (yours is: 1+age+%s), for example: logit(pij)=log(pij/pii)= aij+ bij age + V1 age + etc. as a function of age: <a href=\"%s_%d-1.svg\">%s_%d-1.svg</a><br> \ |
| <img src=\"%s_%d-1.svg\">",model,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1); |
<img src=\"%s_%d-1.svg\">",model,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1); |
| /* Pij */ |
/* Pij */ |
| fprintf(fichtm,"<br>\n- P<sub>ij</sub> or conditional probabilities to be observed in state j being in state i, %d (stepm) months before: <a href=\"%s_%d-2.svg\">%s_%d-2.svg</a><br> \ |
fprintf(fichtm,"<br>\n- P<sub>ij</sub> or conditional probabilities to be observed in state j being in state i, %d (stepm) months before: <a href=\"%s_%d-2.svg\">%s_%d-2.svg</a><br> \ |
| <img src=\"%s_%d-2.svg\">",stepm,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1); |
<img src=\"%s_%d-2.svg\">",stepm,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1); |
| /* Quasi-incidences */ |
/* Quasi-incidences */ |
| fprintf(fichtm,"<br>\n- I<sub>ij</sub> or Conditional probabilities to be observed in state j being in state i %d (stepm) months\ |
fprintf(fichtm,"<br>\n- I<sub>ij</sub> or Conditional probabilities to be observed in state j being in state i %d (stepm) months\ |
| before but expressed in per year i.e. quasi incidences if stepm is small and probabilities too, \ |
before but expressed in per year i.e. quasi incidences if stepm is small and probabilities too, \ |
| incidence (rates) are the limit when h tends to zero of the ratio of the probability <sub>h</sub>P<sub>ij</sub> \ |
incidence (rates) are the limit when h tends to zero of the ratio of the probability <sub>h</sub>P<sub>ij</sub> \ |
| divided by h: <sub>h</sub>P<sub>ij</sub>/h : <a href=\"%s_%d-3.svg\">%s_%d-3.svg</a><br> \ |
divided by h: <sub>h</sub>P<sub>ij</sub>/h : <a href=\"%s_%d-3.svg\">%s_%d-3.svg</a><br> \ |
| <img src=\"%s_%d-3.svg\">",stepm,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1); |
<img src=\"%s_%d-3.svg\">",stepm,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1,subdirf2(optionfilefiname,"PE_"),jj1); |
| /* Survival functions (period) in state j */ |
/* Survival functions (period) in state j */ |
| for(cpt=1; cpt<=nlstate;cpt++){ |
for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Survival functions in state %d. Or probability to survive in state %d being in state (1 to %d) at different ages. <a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> \ |
fprintf(fichtm,"<br>\n- Survival functions in state %d. Or probability to survive in state %d being in state (1 to %d) at different ages. <a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> \ |
| <img src=\"%s_%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"LIJ_"),cpt,jj1,subdirf2(optionfilefiname,"LIJ_"),cpt,jj1,subdirf2(optionfilefiname,"LIJ_"),cpt,jj1); |
<img src=\"%s_%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"LIJ_"),cpt,jj1,subdirf2(optionfilefiname,"LIJ_"),cpt,jj1,subdirf2(optionfilefiname,"LIJ_"),cpt,jj1); |
| } |
} |
| /* State specific survival functions (period) */ |
/* State specific survival functions (period) */ |
| for(cpt=1; cpt<=nlstate;cpt++){ |
for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Survival functions from state %d in each live state and total.\ |
fprintf(fichtm,"<br>\n- Survival functions from state %d in each live state and total.\ |
| Or probability to survive in various states (1 to %d) being in state %d at different ages. \ |
Or probability to survive in various states (1 to %d) being in state %d at different ages. \ |
| <a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> <img src=\"%s_%d-%d.svg\">", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,jj1,subdirf2(optionfilefiname,"LIJT_"),cpt,jj1,subdirf2(optionfilefiname,"LIJT_"),cpt,jj1); |
<a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> <img src=\"%s_%d-%d.svg\">", cpt, nlstate, cpt, subdirf2(optionfilefiname,"LIJT_"),cpt,jj1,subdirf2(optionfilefiname,"LIJT_"),cpt,jj1,subdirf2(optionfilefiname,"LIJT_"),cpt,jj1); |
| } |
} |
| /* Period (stable) prevalence in each health state */ |
/* Period (stable) prevalence in each health state */ |
| for(cpt=1; cpt<=nlstate;cpt++){ |
|
| fprintf(fichtm,"<br>\n- Convergence to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s_%d-%d.svg\">%s_%d-%d.svg</a><br> \ |
|
| <img src=\"%s_%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"P_"),cpt,jj1,subdirf2(optionfilefiname,"P_"),cpt,jj1,subdirf2(optionfilefiname,"P_"),cpt,jj1); |
|
| } |
|
| if(backcast==1){ |
|
| /* Period (stable) back prevalence in each health state */ |
|
| for(cpt=1; cpt<=nlstate;cpt++){ |
for(cpt=1; cpt<=nlstate;cpt++){ |
| fprintf(fichtm,"<br>\n- Convergence to period (stable) back prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s_%d-%d.svg\">%s_%d-%d.svg</a><br> \ |
fprintf(fichtm,"<br>\n- Convergence to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s_%d-%d.svg\">%s_%d-%d.svg</a><br> \ |
| |
<img src=\"%s_%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"P_"),cpt,jj1,subdirf2(optionfilefiname,"P_"),cpt,jj1,subdirf2(optionfilefiname,"P_"),cpt,jj1); |
| |
} |
| |
if(backcast==1){ |
| |
/* Period (stable) back prevalence in each health state */ |
| |
for(cpt=1; cpt<=nlstate;cpt++){ |
| |
fprintf(fichtm,"<br>\n- Convergence to period (stable) back prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s_%d-%d.svg\">%s_%d-%d.svg</a><br> \ |
| <img src=\"%s_%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"PB_"),cpt,jj1,subdirf2(optionfilefiname,"PB_"),cpt,jj1,subdirf2(optionfilefiname,"PB_"),cpt,jj1); |
<img src=\"%s_%d-%d.svg\">", cpt, cpt, nlstate, subdirf2(optionfilefiname,"PB_"),cpt,jj1,subdirf2(optionfilefiname,"PB_"),cpt,jj1,subdirf2(optionfilefiname,"PB_"),cpt,jj1); |
| |
} |
| } |
} |
| } |
if(prevfcast==1){ |
| if(prevfcast==1){ |
/* Projection of prevalence up to period (stable) prevalence in each health state */ |
| /* Projection of prevalence up to period (stable) prevalence in each health state */ |
for(cpt=1; cpt<=nlstate;cpt++){ |
| for(cpt=1; cpt<=nlstate;cpt++){ |
fprintf(fichtm,"<br>\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f) up to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> \ |
| fprintf(fichtm,"<br>\n- Projection of cross-sectional prevalence (estimated with cases observed from %.1f to %.1f) up to period (stable) prevalence in state %d. Or probability to be in state %d being in state (1 to %d) at different ages. <a href=\"%s%d_%d.svg\">%s%d_%d.svg</a><br> \ |
|
| <img src=\"%s_%d-%d.svg\">", dateprev1, dateprev2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1); |
<img src=\"%s_%d-%d.svg\">", dateprev1, dateprev2, cpt, cpt, nlstate, subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1,subdirf2(optionfilefiname,"PROJ_"),cpt,jj1); |
| } |
} |
| } |
} |
| |
|
| for(cpt=1; cpt<=nlstate;cpt++) { |
for(cpt=1; cpt<=nlstate;cpt++) { |
| fprintf(fichtm,"\n<br>- Life expectancy by health state (%d) at initial age and its decomposition into health expectancies in each alive state (1 to %d) (or area under each survival functions): <a href=\"%s_%d%d.svg\">%s_%d%d.svg</a> <br> \ |
fprintf(fichtm,"\n<br>- Life expectancy by health state (%d) at initial age and its decomposition into health expectancies in each alive state (1 to %d) (or area under each survival functions): <a href=\"%s_%d%d.svg\">%s_%d%d.svg</a> <br> \ |
| <img src=\"%s_%d%d.svg\">",cpt,nlstate,subdirf2(optionfilefiname,"EXP_"),cpt,jj1,subdirf2(optionfilefiname,"EXP_"),cpt,jj1,subdirf2(optionfilefiname,"EXP_"),cpt,jj1); |
<img src=\"%s_%d%d.svg\">",cpt,nlstate,subdirf2(optionfilefiname,"EXP_"),cpt,jj1,subdirf2(optionfilefiname,"EXP_"),cpt,jj1,subdirf2(optionfilefiname,"EXP_"),cpt,jj1); |
| } |
} |
| /* } /\* end i1 *\/ */ |
/* } /\* end i1 *\/ */ |
| }/* End k1 */ |
}/* End k1 */ |
| fprintf(fichtm,"</ul>"); |
fprintf(fichtm,"</ul>"); |
| |
|
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| \n<br><li><h4> <a name='secondorder'>Result files (second order: variances)</a></h4>\n\ |
\n<br><li><h4> <a name='secondorder'>Result files (second order: variances)</a></h4>\n\ |
| - Parameter file with estimated parameters and covariance matrix: <a href=\"%s\">%s</a> <br> \ |
- Parameter file with estimated parameters and covariance matrix: <a href=\"%s\">%s</a> <br> \ |
| - 95%% confidence intervals and Wald tests of the estimated parameters are in the log file if optimization has been done (mle != 0).<br> \ |
- 95%% confidence intervals and Wald tests of the estimated parameters are in the log file if optimization has been done (mle != 0).<br> \ |
|
Line 5622 variances but at the covariance matrix.
|
Line 5625 variances but at the covariance matrix.
|
| covariance matrix of the one-step probabilities. \ |
covariance matrix of the one-step probabilities. \ |
| See page 'Matrix of variance-covariance of one-step probabilities' below. \n", rfileres,rfileres); |
See page 'Matrix of variance-covariance of one-step probabilities' below. \n", rfileres,rfileres); |
| |
|
| fprintf(fichtm," - Standard deviation of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
fprintf(fichtm," - Standard deviation of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
| subdirf2(fileresu,"PROB_"),subdirf2(fileresu,"PROB_")); |
subdirf2(fileresu,"PROB_"),subdirf2(fileresu,"PROB_")); |
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| - Variance-covariance of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
- Variance-covariance of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
| subdirf2(fileresu,"PROBCOV_"),subdirf2(fileresu,"PROBCOV_")); |
subdirf2(fileresu,"PROBCOV_"),subdirf2(fileresu,"PROBCOV_")); |
| |
|
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| - Correlation matrix of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
- Correlation matrix of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
| subdirf2(fileresu,"PROBCOR_"),subdirf2(fileresu,"PROBCOR_")); |
subdirf2(fileresu,"PROBCOR_"),subdirf2(fileresu,"PROBCOR_")); |
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| - Variances and covariances of health expectancies by age and <b>initial health status</b> (cov(e<sup>ij</sup>,e<sup>kl</sup>)(estepm=%2d months): \ |
- Variances and covariances of health expectancies by age and <b>initial health status</b> (cov(e<sup>ij</sup>,e<sup>kl</sup>)(estepm=%2d months): \ |
| <a href=\"%s\">%s</a> <br>\n</li>", |
<a href=\"%s\">%s</a> <br>\n</li>", |
| estepm,subdirf2(fileresu,"CVE_"),subdirf2(fileresu,"CVE_")); |
estepm,subdirf2(fileresu,"CVE_"),subdirf2(fileresu,"CVE_")); |
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| - (a) Health expectancies by health status at initial age (e<sup>ij</sup>) and standard errors (in parentheses) (b) life expectancies and standard errors (e<sup>i.</sup>=e<sup>i1</sup>+e<sup>i2</sup>+...)(estepm=%2d months): \ |
- (a) Health expectancies by health status at initial age (e<sup>ij</sup>) and standard errors (in parentheses) (b) life expectancies and standard errors (e<sup>i.</sup>=e<sup>i1</sup>+e<sup>i2</sup>+...)(estepm=%2d months): \ |
| <a href=\"%s\">%s</a> <br>\n</li>", |
<a href=\"%s\">%s</a> <br>\n</li>", |
| estepm,subdirf2(fileresu,"STDE_"),subdirf2(fileresu,"STDE_")); |
estepm,subdirf2(fileresu,"STDE_"),subdirf2(fileresu,"STDE_")); |
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| - Variances and covariances of health expectancies by age. Status (i) based health expectancies (in state j), e<sup>ij</sup> are weighted by the period prevalences in each state i (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a><br>\n", |
- Variances and covariances of health expectancies by age. Status (i) based health expectancies (in state j), e<sup>ij</sup> are weighted by the period prevalences in each state i (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a><br>\n", |
| estepm, subdirf2(fileresu,"V_"),subdirf2(fileresu,"V_")); |
estepm, subdirf2(fileresu,"V_"),subdirf2(fileresu,"V_")); |
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| - Total life expectancy and total health expectancies to be spent in each health state e<sup>.j</sup> with their standard errors (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a> <br>\n", |
- Total life expectancy and total health expectancies to be spent in each health state e<sup>.j</sup> with their standard errors (if popbased=1, an additional computation is done using the cross-sectional prevalences, i.e population based) (estepm=%d months): <a href=\"%s\">%s</a> <br>\n", |
| estepm, subdirf2(fileresu,"T_"),subdirf2(fileresu,"T_")); |
estepm, subdirf2(fileresu,"T_"),subdirf2(fileresu,"T_")); |
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| - Standard deviation of period (stable) prevalences: <a href=\"%s\">%s</a> <br>\n",\ |
- Standard deviation of period (stable) prevalences: <a href=\"%s\">%s</a> <br>\n",\ |
| subdirf2(fileresu,"VPL_"),subdirf2(fileresu,"VPL_")); |
subdirf2(fileresu,"VPL_"),subdirf2(fileresu,"VPL_")); |
| |
|
| /* if(popforecast==1) fprintf(fichtm,"\n */ |
/* if(popforecast==1) fprintf(fichtm,"\n */ |
| /* - Prevalences forecasting: <a href=\"f%s\">f%s</a> <br>\n */ |
/* - Prevalences forecasting: <a href=\"f%s\">f%s</a> <br>\n */ |
|
Line 5655 See page 'Matrix of variance-covariance
|
Line 5658 See page 'Matrix of variance-covariance
|
| /* <br>",fileres,fileres,fileres,fileres); */ |
/* <br>",fileres,fileres,fileres,fileres); */ |
| /* else */ |
/* else */ |
| /* fprintf(fichtm,"\n No population forecast: popforecast = %d (instead of 1) or stepm = %d (instead of 1) or model=%s (instead of .)<br><br></li>\n",popforecast, stepm, model); */ |
/* fprintf(fichtm,"\n No population forecast: popforecast = %d (instead of 1) or stepm = %d (instead of 1) or model=%s (instead of .)<br><br></li>\n",popforecast, stepm, model); */ |
| fflush(fichtm); |
fflush(fichtm); |
| fprintf(fichtm," <ul><li><b>Graphs</b></li><p>"); |
fprintf(fichtm," <ul><li><b>Graphs</b></li><p>"); |
| |
|
| m=pow(2,cptcoveff); |
m=pow(2,cptcoveff); |
| if (cptcovn < 1) {m=1;ncodemax[1]=1;} |
if (cptcovn < 1) {m=1;ncodemax[1]=1;} |
| |
|
| jj1=0; |
jj1=0; |
| for(k1=1; k1<=m;k1++){ |
for(k1=1; k1<=m;k1++){ |
| /* for(i1=1; i1<=ncodemax[k1];i1++){ */ |
/* for(i1=1; i1<=ncodemax[k1];i1++){ */ |
| jj1++; |
jj1++; |
| if (cptcovn > 0) { |
if (cptcovn > 0) { |
| fprintf(fichtm,"<hr size=\"2\" color=\"#EC5E5E\">************ Results for covariates"); |
fprintf(fichtm,"<hr size=\"2\" color=\"#EC5E5E\">************ Results for covariates"); |
| for (cpt=1; cpt<=cptcoveff;cpt++) |
for (cpt=1; cpt<=cptcoveff;cpt++) |
| fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]); |
fprintf(fichtm," V%d=%d ",Tvaraff[cpt],nbcode[Tvaraff[cpt]][codtabm(jj1,cpt)]); |
| fprintf(fichtm," ************\n<hr size=\"2\" color=\"#EC5E5E\">"); |
fprintf(fichtm," ************\n<hr size=\"2\" color=\"#EC5E5E\">"); |
| |
|
| if(invalidvarcomb[k1]){ |
if(invalidvarcomb[k1]){ |
| fprintf(fichtm,"\n<h4>Combination (%d) ignored because no cases </h4>\n",k1); |
fprintf(fichtm,"\n<h4>Combination (%d) ignored because no cases </h4>\n",k1); |
| continue; |
continue; |
| } |
} |
| } |
} |
| for(cpt=1; cpt<=nlstate;cpt++) { |
for(cpt=1; cpt<=nlstate;cpt++) { |
| fprintf(fichtm,"\n<br>- Observed (cross-sectional) and period (incidence based) \ |
fprintf(fichtm,"\n<br>- Observed (cross-sectional) and period (incidence based) \ |
|
Line 5687 true period expectancies (those weighted
|
Line 5690 true period expectancies (those weighted
|
| drawn in addition to the population based expectancies computed using\ |
drawn in addition to the population based expectancies computed using\ |
| observed and cahotic prevalences: <a href=\"%s_%d.svg\">%s_%d.svg</a>\n<br>\ |
observed and cahotic prevalences: <a href=\"%s_%d.svg\">%s_%d.svg</a>\n<br>\ |
| <img src=\"%s_%d.svg\">",subdirf2(optionfilefiname,"E_"),jj1,subdirf2(optionfilefiname,"E_"),jj1,subdirf2(optionfilefiname,"E_"),jj1); |
<img src=\"%s_%d.svg\">",subdirf2(optionfilefiname,"E_"),jj1,subdirf2(optionfilefiname,"E_"),jj1,subdirf2(optionfilefiname,"E_"),jj1); |
| /* } /\* end i1 *\/ */ |
/* } /\* end i1 *\/ */ |
| }/* End k1 */ |
}/* End k1 */ |
| fprintf(fichtm,"</ul>"); |
fprintf(fichtm,"</ul>"); |
| fflush(fichtm); |
fflush(fichtm); |
| } |
} |
| |
|
| /******************* Gnuplot file **************/ |
/******************* Gnuplot file **************/ |
|
Line 6324 plot [%.f:%.f] ", ageminpar, agemaxpar)
|
Line 6327 plot [%.f:%.f] ", ageminpar, agemaxpar)
|
| |
|
| /*************** Moving average **************/ |
/*************** Moving average **************/ |
| /* int movingaverage(double ***probs, double bage, double fage, double ***mobaverage, int mobilav, double bageout, double fageout){ */ |
/* int movingaverage(double ***probs, double bage, double fage, double ***mobaverage, int mobilav, double bageout, double fageout){ */ |
| int movingaverage(double ***probs, double bage, double fage, double ***mobaverage, int mobilav){ |
int movingaverage(double ***probs, double bage, double fage, double ***mobaverage, int mobilav){ |
| |
|
| int i, cpt, cptcod; |
int i, cpt, cptcod; |
| int modcovmax =1; |
int modcovmax =1; |
| int mobilavrange, mob; |
int mobilavrange, mob; |
| int iage=0; |
int iage=0; |
| |
|
| double sum=0.; |
double sum=0.; |
| double age; |
double age; |
| double *sumnewp, *sumnewm; |
double *sumnewp, *sumnewm; |
| double *agemingood, *agemaxgood; /* Currently identical for all covariates */ |
double *agemingood, *agemaxgood; /* Currently identical for all covariates */ |
| |
|
| |
|
| /* modcovmax=2*cptcoveff;/\* Max number of modalities. We suppose */ |
/* modcovmax=2*cptcoveff;/\* Max number of modalities. We suppose */ |
| /* a covariate has 2 modalities, should be equal to ncovcombmax *\/ */ |
/* a covariate has 2 modalities, should be equal to ncovcombmax *\/ */ |
| |
|
| sumnewp = vector(1,ncovcombmax); |
sumnewp = vector(1,ncovcombmax); |
| sumnewm = vector(1,ncovcombmax); |
sumnewm = vector(1,ncovcombmax); |
| agemingood = vector(1,ncovcombmax); |
agemingood = vector(1,ncovcombmax); |
| agemaxgood = vector(1,ncovcombmax); |
agemaxgood = vector(1,ncovcombmax); |
| |
|
| for (cptcod=1;cptcod<=ncovcombmax;cptcod++){ |
for (cptcod=1;cptcod<=ncovcombmax;cptcod++){ |
| sumnewm[cptcod]=0.; |
sumnewm[cptcod]=0.; |
| sumnewp[cptcod]=0.; |
sumnewp[cptcod]=0.; |
| agemingood[cptcod]=0; |
agemingood[cptcod]=0; |
| agemaxgood[cptcod]=0; |
agemaxgood[cptcod]=0; |
| } |
} |
| if (cptcovn<1) ncovcombmax=1; /* At least 1 pass */ |
if (cptcovn<1) ncovcombmax=1; /* At least 1 pass */ |
| |
|
| if(mobilav==1||mobilav ==3 ||mobilav==5 ||mobilav== 7){ |
if(mobilav==1||mobilav ==3 ||mobilav==5 ||mobilav== 7){ |
| if(mobilav==1) mobilavrange=5; /* default */ |
if(mobilav==1) mobilavrange=5; /* default */ |
| else mobilavrange=mobilav; |
else mobilavrange=mobilav; |
| for (age=bage; age<=fage; age++) |
for (age=bage; age<=fage; age++) |
| for (i=1; i<=nlstate;i++) |
for (i=1; i<=nlstate;i++) |
| for (cptcod=1;cptcod<=ncovcombmax;cptcod++) |
for (cptcod=1;cptcod<=ncovcombmax;cptcod++) |
| mobaverage[(int)age][i][cptcod]=probs[(int)age][i][cptcod]; |
mobaverage[(int)age][i][cptcod]=probs[(int)age][i][cptcod]; |
| /* We keep the original values on the extreme ages bage, fage and for |
/* We keep the original values on the extreme ages bage, fage and for |
| fage+1 and bage-1 we use a 3 terms moving average; for fage+2 bage+2 |
fage+1 and bage-1 we use a 3 terms moving average; for fage+2 bage+2 |
| we use a 5 terms etc. until the borders are no more concerned. |
we use a 5 terms etc. until the borders are no more concerned. |
| */ |
*/ |
| for (mob=3;mob <=mobilavrange;mob=mob+2){ |
for (mob=3;mob <=mobilavrange;mob=mob+2){ |
| for (age=bage+(mob-1)/2; age<=fage-(mob-1)/2; age++){ |
for (age=bage+(mob-1)/2; age<=fage-(mob-1)/2; age++){ |
| for (i=1; i<=nlstate;i++){ |
for (i=1; i<=nlstate;i++){ |
| for (cptcod=1;cptcod<=ncovcombmax;cptcod++){ |
for (cptcod=1;cptcod<=ncovcombmax;cptcod++){ |
| mobaverage[(int)age][i][cptcod] =probs[(int)age][i][cptcod]; |
mobaverage[(int)age][i][cptcod] =probs[(int)age][i][cptcod]; |
| for (cpt=1;cpt<=(mob-1)/2;cpt++){ |
for (cpt=1;cpt<=(mob-1)/2;cpt++){ |
| mobaverage[(int)age][i][cptcod] +=probs[(int)age-cpt][i][cptcod]; |
mobaverage[(int)age][i][cptcod] +=probs[(int)age-cpt][i][cptcod]; |
| mobaverage[(int)age][i][cptcod] +=probs[(int)age+cpt][i][cptcod]; |
mobaverage[(int)age][i][cptcod] +=probs[(int)age+cpt][i][cptcod]; |
| } |
} |
| mobaverage[(int)age][i][cptcod]=mobaverage[(int)age][i][cptcod]/mob; |
mobaverage[(int)age][i][cptcod]=mobaverage[(int)age][i][cptcod]/mob; |
| } |
} |
| } |
} |
| }/* end age */ |
}/* end age */ |
| }/* end mob */ |
}/* end mob */ |
| }else |
}else |
| return -1; |
return -1; |
| for (cptcod=1;cptcod<=ncovcombmax;cptcod++){ |
for (cptcod=1;cptcod<=ncovcombmax;cptcod++){ |
| /* for (age=bage+(mob-1)/2; age<=fage-(mob-1)/2; age++){ */ |
/* for (age=bage+(mob-1)/2; age<=fage-(mob-1)/2; age++){ */ |
| agemingood[cptcod]=fage-(mob-1)/2; |
if(invalidvarcomb[cptcod]){ |
| for (age=fage-(mob-1)/2; age>=bage; age--){/* From oldest to youngest, finding the youngest wrong */ |
printf("\nCombination (%d) ignored because no cases \n",cptcod); |
| sumnewm[cptcod]=0.; |
continue; |
| for (i=1; i<=nlstate;i++){ |
} |
| sumnewm[cptcod]+=mobaverage[(int)age][i][cptcod]; |
|
| } |
agemingood[cptcod]=fage-(mob-1)/2; |
| if(fabs(sumnewm[cptcod] - 1.) <= 1.e-3) { /* good */ |
for (age=fage-(mob-1)/2; age>=bage; age--){/* From oldest to youngest, finding the youngest wrong */ |
| agemingood[cptcod]=age; |
sumnewm[cptcod]=0.; |
| }else{ /* bad */ |
for (i=1; i<=nlstate;i++){ |
| for (i=1; i<=nlstate;i++){ |
sumnewm[cptcod]+=mobaverage[(int)age][i][cptcod]; |
| mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; |
} |
| } /* i */ |
if(fabs(sumnewm[cptcod] - 1.) <= 1.e-3) { /* good */ |
| } /* end bad */ |
agemingood[cptcod]=age; |
| }/* age */ |
}else{ /* bad */ |
| sum=0.; |
for (i=1; i<=nlstate;i++){ |
| for (i=1; i<=nlstate;i++){ |
mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; |
| sum+=mobaverage[(int)agemingood[cptcod]][i][cptcod]; |
} /* i */ |
| } |
} /* end bad */ |
| if(fabs(sum - 1.) > 1.e-3) { /* bad */ |
}/* age */ |
| printf("For this combination of covariate cptcod=%d, we can't get a smoothed prevalence which sums to one at any descending age!\n",cptcod); |
sum=0.; |
| /* for (i=1; i<=nlstate;i++){ */ |
for (i=1; i<=nlstate;i++){ |
| /* mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; */ |
sum+=mobaverage[(int)agemingood[cptcod]][i][cptcod]; |
| /* } /\* i *\/ */ |
} |
| } /* end bad */ |
if(fabs(sum - 1.) > 1.e-3) { /* bad */ |
| /* else{ /\* We found some ages summing to one, we will smooth the oldest *\/ */ |
printf("For this combination of covariate cptcod=%d, we can't get a smoothed prevalence which sums to one at any descending age!\n",cptcod); |
| /* From youngest, finding the oldest wrong */ |
/* for (i=1; i<=nlstate;i++){ */ |
| agemaxgood[cptcod]=bage+(mob-1)/2; |
/* mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; */ |
| for (age=bage+(mob-1)/2; age<=fage; age++){ |
/* } /\* i *\/ */ |
| sumnewm[cptcod]=0.; |
} /* end bad */ |
| for (i=1; i<=nlstate;i++){ |
/* else{ /\* We found some ages summing to one, we will smooth the oldest *\/ */ |
| sumnewm[cptcod]+=mobaverage[(int)age][i][cptcod]; |
/* From youngest, finding the oldest wrong */ |
| } |
agemaxgood[cptcod]=bage+(mob-1)/2; |
| if(fabs(sumnewm[cptcod] - 1.) <= 1.e-3) { /* good */ |
for (age=bage+(mob-1)/2; age<=fage; age++){ |
| agemaxgood[cptcod]=age; |
sumnewm[cptcod]=0.; |
| }else{ /* bad */ |
for (i=1; i<=nlstate;i++){ |
| for (i=1; i<=nlstate;i++){ |
sumnewm[cptcod]+=mobaverage[(int)age][i][cptcod]; |
| mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemaxgood[cptcod]][i][cptcod]; |
} |
| } /* i */ |
if(fabs(sumnewm[cptcod] - 1.) <= 1.e-3) { /* good */ |
| } /* end bad */ |
agemaxgood[cptcod]=age; |
| }/* age */ |
}else{ /* bad */ |
| sum=0.; |
for (i=1; i<=nlstate;i++){ |
| for (i=1; i<=nlstate;i++){ |
mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemaxgood[cptcod]][i][cptcod]; |
| sum+=mobaverage[(int)agemaxgood[cptcod]][i][cptcod]; |
} /* i */ |
| } |
} /* end bad */ |
| if(fabs(sum - 1.) > 1.e-3) { /* bad */ |
}/* age */ |
| printf("For this combination of covariate cptcod=%d, we can't get a smoothed prevalence which sums to one at any ascending age!\n",cptcod); |
sum=0.; |
| /* for (i=1; i<=nlstate;i++){ */ |
for (i=1; i<=nlstate;i++){ |
| /* mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; */ |
sum+=mobaverage[(int)agemaxgood[cptcod]][i][cptcod]; |
| /* } /\* i *\/ */ |
} |
| } /* end bad */ |
if(fabs(sum - 1.) > 1.e-3) { /* bad */ |
| |
printf("For this combination of covariate cptcod=%d, we can't get a smoothed prevalence which sums to one at any ascending age!\n",cptcod); |
| for (age=bage; age<=fage; age++){ |
/* for (i=1; i<=nlstate;i++){ */ |
| printf("%d %d ", cptcod, (int)age); |
/* mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; */ |
| sumnewp[cptcod]=0.; |
/* } /\* i *\/ */ |
| sumnewm[cptcod]=0.; |
} /* end bad */ |
| for (i=1; i<=nlstate;i++){ |
|
| sumnewp[cptcod]+=probs[(int)age][i][cptcod]; |
for (age=bage; age<=fage; age++){ |
| sumnewm[cptcod]+=mobaverage[(int)age][i][cptcod]; |
printf("%d %d ", cptcod, (int)age); |
| /* printf("%.4f %.4f ",probs[(int)age][i][cptcod], mobaverage[(int)age][i][cptcod]); */ |
sumnewp[cptcod]=0.; |
| } |
sumnewm[cptcod]=0.; |
| /* printf("%.4f %.4f \n",sumnewp[cptcod], sumnewm[cptcod]); */ |
for (i=1; i<=nlstate;i++){ |
| } |
sumnewp[cptcod]+=probs[(int)age][i][cptcod]; |
| /* printf("\n"); */ |
sumnewm[cptcod]+=mobaverage[(int)age][i][cptcod]; |
| /* } */ |
/* printf("%.4f %.4f ",probs[(int)age][i][cptcod], mobaverage[(int)age][i][cptcod]); */ |
| /* brutal averaging */ |
} |
| for (i=1; i<=nlstate;i++){ |
/* printf("%.4f %.4f \n",sumnewp[cptcod], sumnewm[cptcod]); */ |
| for (age=1; age<=bage; age++){ |
} |
| mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; |
/* printf("\n"); */ |
| /* printf("age=%d i=%d cptcod=%d mobaverage=%.4f \n",(int)age,i, cptcod, mobaverage[(int)age][i][cptcod]); */ |
/* } */ |
| } |
/* brutal averaging */ |
| for (age=fage; age<=AGESUP; age++){ |
for (i=1; i<=nlstate;i++){ |
| mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemaxgood[cptcod]][i][cptcod]; |
for (age=1; age<=bage; age++){ |
| /* printf("age=%d i=%d cptcod=%d mobaverage=%.4f \n",(int)age,i, cptcod, mobaverage[(int)age][i][cptcod]); */ |
mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemingood[cptcod]][i][cptcod]; |
| } |
/* printf("age=%d i=%d cptcod=%d mobaverage=%.4f \n",(int)age,i, cptcod, mobaverage[(int)age][i][cptcod]); */ |
| } /* end i status */ |
} |
| for (i=nlstate+1; i<=nlstate+ndeath;i++){ |
for (age=fage; age<=AGESUP; age++){ |
| for (age=1; age<=AGESUP; age++){ |
mobaverage[(int)age][i][cptcod]=mobaverage[(int)agemaxgood[cptcod]][i][cptcod]; |
| /*printf("i=%d, age=%d, cptcod=%d\n",i, (int)age, cptcod);*/ |
/* printf("age=%d i=%d cptcod=%d mobaverage=%.4f \n",(int)age,i, cptcod, mobaverage[(int)age][i][cptcod]); */ |
| mobaverage[(int)age][i][cptcod]=0.; |
} |
| } |
} /* end i status */ |
| } |
for (i=nlstate+1; i<=nlstate+ndeath;i++){ |
| }/* end cptcod */ |
for (age=1; age<=AGESUP; age++){ |
| free_vector(sumnewm,1, ncovcombmax); |
/*printf("i=%d, age=%d, cptcod=%d\n",i, (int)age, cptcod);*/ |
| free_vector(sumnewp,1, ncovcombmax); |
mobaverage[(int)age][i][cptcod]=0.; |
| free_vector(agemaxgood,1, ncovcombmax); |
} |
| free_vector(agemingood,1, ncovcombmax); |
} |
| return 0; |
}/* end cptcod */ |
| }/* End movingaverage */ |
free_vector(sumnewm,1, ncovcombmax); |
| |
free_vector(sumnewp,1, ncovcombmax); |
| |
free_vector(agemaxgood,1, ncovcombmax); |
| |
free_vector(agemingood,1, ncovcombmax); |
| |
return 0; |
| |
}/* End movingaverage */ |
| |
|
| |
|
| /************** Forecasting ******************/ |
/************** Forecasting ******************/ |
|
Line 8131 int hPijx(double *p, int bage, int fage)
|
Line 8139 int hPijx(double *p, int bage, int fage)
|
| for(j=1;j<=cptcoveff;j++) |
for(j=1;j<=cptcoveff;j++) |
| fprintf(ficrespijb,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]); |
fprintf(ficrespijb,"V%d=%d ",Tvaraff[j],nbcode[Tvaraff[j]][codtabm(k,j)]); |
| fprintf(ficrespijb,"******\n"); |
fprintf(ficrespijb,"******\n"); |
| |
if(invalidvarcomb[k]){ |
| |
fprintf(ficrespijb,"\n#Combination (%d) ignored because no cases \n",k); |
| |
continue; |
| |
} |
| |
|
| /* for (agedeb=fage; agedeb>=bage; agedeb--){ /\* If stepm=6 months *\/ */ |
/* for (agedeb=fage; agedeb>=bage; agedeb--){ /\* If stepm=6 months *\/ */ |
| for (agedeb=bage; agedeb<=fage; agedeb++){ /* If stepm=6 months and estepm=24 (2 years) */ |
for (agedeb=bage; agedeb<=fage; agedeb++){ /* If stepm=6 months and estepm=24 (2 years) */ |
|
Line 9648 Please run with mle=-1 to get a correct
|
Line 9660 Please run with mle=-1 to get a correct
|
| back_prevalence_limit(p, bprlim, ageminpar, agemaxpar, ftolpl, &ncvyear, dateprev1, dateprev2, firstpass, lastpass, mobilavproj); |
back_prevalence_limit(p, bprlim, ageminpar, agemaxpar, ftolpl, &ncvyear, dateprev1, dateprev2, firstpass, lastpass, mobilavproj); |
| fclose(ficresplb); |
fclose(ficresplb); |
| |
|
| /* hBijx(p, bage, fage, mobaverage); */ |
hBijx(p, bage, fage, mobaverage); |
| /* fclose(ficrespijb); */ |
fclose(ficrespijb); |
| free_matrix(bprlim,1,nlstate,1,nlstate); /*here or after loop ? */ |
free_matrix(bprlim,1,nlstate,1,nlstate); /*here or after loop ? */ |
| |
|
| /* prevbackforecast(fileresu, anback1, mback1, jback1, agemin, agemax, dateprev1, dateprev2, mobilavproj, |
/* prevbackforecast(fileresu, anback1, mback1, jback1, agemin, agemax, dateprev1, dateprev2, mobilavproj, |
![[BACK]](/icons/small/back.gif){kind=icon}
Return to imach.c CVS log ![[TXT]](/icons/small/text.gif){kind=icon}
![[DIR]](/icons/small/dir.gif){kind=icon}
Up to [Local Repository] / imach / src