version 1.50, 2002/06/26 23:25:02
|
version 1.52, 2002/07/19 18:49:30
|
Line 1862 void varevsij(char optionfilefiname[], d
|
Line 1862 void varevsij(char optionfilefiname[], d
|
exit(0);
|
exit(0);
|
}
|
}
|
else{
|
else{
|
fprintf(fichtm,"\n<li><h4> Computing step probabilities of dying and weighted average (i.e global mortality independent of initial healh state)</h4></li>\n");
|
fprintf(fichtm,"\n<li><h4> Computing probabilities of dying as a weighted average (i.e global mortality independent of initial healh state)</h4></li>\n");
|
}
|
}
|
varppt = matrix(nlstate+1,nlstate+ndeath,nlstate+1,nlstate+ndeath);
|
varppt = matrix(nlstate+1,nlstate+ndeath,nlstate+1,nlstate+ndeath);
|
|
|
Line 2013 void varevsij(char optionfilefiname[], d
|
Line 2013 void varevsij(char optionfilefiname[], d
|
/* This for computing force of mortality (h=1)as a weighted average */
|
/* This for computing force of mortality (h=1)as a weighted average */
|
for(j=nlstate+1,gmp[j]=0.;j<=nlstate+ndeath;j++){
|
for(j=nlstate+1,gmp[j]=0.;j<=nlstate+ndeath;j++){
|
for(i=1; i<= nlstate; i++)
|
for(i=1; i<= nlstate; i++)
|
gmp[j] += prlim[i][i]*p3mat[i][j][1];
|
gmp[j] += prlim[i][i]*p3mat[i][j][1];
|
}
|
}
|
/* end force of mortality */
|
/* end force of mortality */
|
|
|
Line 2049 void varevsij(char optionfilefiname[], d
|
Line 2049 void varevsij(char optionfilefiname[], d
|
fprintf(ficgp,"\n replot \"%s\" u 1:(($3+1.96*$4)*%6.3f) t \"95\%% interval\" w l 2 ",fileresprobmorprev,YEARM/estepm);
|
fprintf(ficgp,"\n replot \"%s\" u 1:(($3+1.96*$4)*%6.3f) t \"95\%% interval\" w l 2 ",fileresprobmorprev,YEARM/estepm);
|
fprintf(ficgp,"\n replot \"%s\" u 1:(($3-1.96*$4)*%6.3f) not w l 2 ",fileresprobmorprev,YEARM/estepm);
|
fprintf(ficgp,"\n replot \"%s\" u 1:(($3-1.96*$4)*%6.3f) not w l 2 ",fileresprobmorprev,YEARM/estepm);
|
fprintf(fichtm,"\n<br> File (multiple files are possible if covariates are present): <A href=\"%s\">%s</a>\n",fileresprobmorprev,fileresprobmorprev);
|
fprintf(fichtm,"\n<br> File (multiple files are possible if covariates are present): <A href=\"%s\">%s</a>\n",fileresprobmorprev,fileresprobmorprev);
|
fprintf(fichtm,"\n<br> Probability is computed over estepm=%d months and then divided by estepm and multiplied by %.0f in order to have the probability to die over a year <br> <img src=\"varmuptjgr%s%s.png\"> <br>\n", stepm,YEARM,digitp,digit);
|
fprintf(fichtm,"\n<br> Probability is computed over estepm=%d months. <br> <img src=\"varmuptjgr%s%s.png\"> <br>\n", stepm,digitp,digit);
|
|
/* fprintf(fichtm,"\n<br> Probability is computed over estepm=%d months and then divided by estepm and multiplied by %.0f in order to have the probability to die over a year <br> <img src=\"varmuptjgr%s%s.png\"> <br>\n", stepm,YEARM,digitp,digit);
|
|
*/
|
fprintf(ficgp,"\nset out \"varmuptjgr%s%s.png\";replot;",digitp,digit);
|
fprintf(ficgp,"\nset out \"varmuptjgr%s%s.png\";replot;",digitp,digit);
|
|
|
free_vector(xp,1,npar);
|
free_vector(xp,1,npar);
|
Line 2154 void varprob(char optionfilefiname[], do
|
Line 2156 void varprob(char optionfilefiname[], do
|
int k2, l2, j1, z1;
|
int k2, l2, j1, z1;
|
int k=0,l, cptcode;
|
int k=0,l, cptcode;
|
int first=1, first1;
|
int first=1, first1;
|
double cv12, mu1, mu2, lc1, lc2, v12, v21, v11, v22,v1,v2;
|
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;
|
Line 2378 void varprob(char optionfilefiname[], do
|
Line 2380 void varprob(char optionfilefiname[], do
|
|
|
/* 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;
|
first1=1;
|
for (k1=1; k1<=(nlstate);k1++){
|
for (k2=1; k2<=(nlstate);k2++){
|
for (l1=1; l1<=(nlstate+ndeath);l1++){
|
for (l2=1; l2<=(nlstate+ndeath);l2++){
|
if(l1==k1) continue;
|
if(l2==k2) continue;
|
i=(k1-1)*(nlstate+ndeath)+l1;
|
j=(k2-1)*(nlstate+ndeath)+l2;
|
for (k2=1; k2<=(nlstate);k2++){
|
for (k1=1; k1<=(nlstate);k1++){
|
for (l2=1; l2<=(nlstate+ndeath);l2++){
|
for (l1=1; l1<=(nlstate+ndeath);l1++){
|
if(l2==k2) continue;
|
if(l1==k1) continue;
|
j=(k2-1)*(nlstate+ndeath)+l2;
|
i=(k1-1)*(nlstate+ndeath)+l1;
|
if(j<=i) 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;
|
Line 2394 void varprob(char optionfilefiname[], do
|
Line 2396 void varprob(char optionfilefiname[], do
|
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);
|
/* 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));
|
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));
|
lc2=((v1+v2)-sqrt((v1+v2)*(v1+v2) - 4*(v1*v2-cv12*cv12)))/2.;
|
if(first1==1){
|
|
first1=0;
|
|
printf("Var %.4e %.4e cov %.4e Eigen %.3e %.3e\nOthers in log...\n",v1,v2,cv12,lc1,lc2);
|
|
}
|
|
fprintf(ficlog,"Var %.4e %.4e cov %.4e Eigen %.3e %.3e\n",v1,v2,cv12,lc1,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);
|
/*v21=sqrt(1.-v11*v11); *//* error */
|
|
v21=(lc1-v1)/cv12*v11;
|
v12=-v21;
|
v12=-v21;
|
v22=v11;
|
v22=v11;
|
|
tnalp=v21/v11;
|
|
if(first1==1){
|
|
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);
|
|
}
|
|
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 + v21*lc2*sin(t) */
|
/* mu2+ v21*lc1*cost + v22*lc2*sin(t) */
|
if(first==1){
|
if(first==1){
|
first=0;
|
first=0;
|
fprintf(ficgp,"\nset parametric;set nolabel");
|
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)\"",k2,l2,k1,l1);
|
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 png small\nset size 0.65,0.65");
|
fprintf(ficgp,"\nset ter png small\nset size 0.65,0.65");
|
fprintf(fichtm,"\n<br>Ellipsoids of confidence cov(p%1d%1d,p%1d%1d) expressed in year<sup>-1</sup> :<a href=\"varpijgr%s%d%1d%1d-%1d%1d.png\">varpijgr%s%d%1d%1d-%1d%1d.png</A>, ",k2,l2,k1,l1,optionfilefiname, j1,k2,l2,k1,l1,optionfilefiname, j1,k2,l2,k1,l1);
|
fprintf(fichtm,"\n<br>Ellipsoids of confidence cov(p%1d%1d,p%1d%1d) expressed in year<sup>-1</sup> :<a href=\"varpijgr%s%d%1d%1d-%1d%1d.png\">varpijgr%s%d%1d%1d-%1d%1d.png</A>, ",k1,l1,k2,l2,optionfilefiname, j1,k1,l1,k2,l2,optionfilefiname, j1,k1,l1,k2,l2);
|
fprintf(fichtm,"\n<br><img src=\"varpijgr%s%d%1d%1d-%1d%1d.png\"> ",optionfilefiname, j1,k2,l2,k1,l1);
|
fprintf(fichtm,"\n<br><img src=\"varpijgr%s%d%1d%1d-%1d%1d.png\"> ",optionfilefiname, j1,k1,l1,k2,l2);
|
fprintf(ficgp,"\nset out \"varpijgr%s%d%1d%1d-%1d%1d.png\"",optionfilefiname, j1,k2,l2,k1,l1);
|
fprintf(fichtm,"\n<br> Correlation at age %d (%.3f),",(int) age, c12);
|
fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu2,mu1);
|
fprintf(ficgp,"\nset out \"varpijgr%s%d%1d%1d-%1d%1d.png\"",optionfilefiname, j1,k1,l1,k2,l2);
|
fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k2,l2,k1,l1);
|
fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu1,mu2);
|
/* 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)) t \"%d\"",\
|
fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k1,l1,k2,l2);
|
mu2,std,v21,sqrt(lc1),v21,sqrt(lc2), \
|
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),(int) age);
|
mu1,std,v11,sqrt(lc1),v12,sqrt(lc2),\
|
*/
|
mu2,std,v21,sqrt(lc1),v22,sqrt(lc2));
|
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",\
|
|
mu2,std,v21,sqrt(lc1),v21,sqrt(lc2), \
|
|
mu1,std,v11,sqrt(lc1),v12,sqrt(lc2));
|
|
}else{
|
}else{
|
first=0;
|
first=0;
|
fprintf(ficgp,"\n# Age %d, p%1d%1d - p%1d%1d",(int) age, k2,l2,k1,l1);
|
fprintf(fichtm," %d (%.3f),",(int) age, c12);
|
fprintf(ficgp,"\nset label \"%d\" at %11.3e,%11.3e center",(int) age, mu2,mu1);
|
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,"\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)) t \"%d\"",\
|
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",\
|
mu2,std,v21,sqrt(lc1),v21,sqrt(lc2), \
|
mu1,std,v11,sqrt(lc1),v12,sqrt(lc2),\
|
mu1,std,v11,sqrt(lc1),v12,sqrt(lc2),(int) age);
|
mu2,std,v21,sqrt(lc1),v22,sqrt(lc2));
|
*/
|
|
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",\
|
|
mu2,std,v21,sqrt(lc1),v21,sqrt(lc2), \
|
|
mu1,std,v11,sqrt(lc1),v12,sqrt(lc2));
|
|
}/* if first */
|
}/* if first */
|
} /* age mod 5 */
|
} /* age mod 5 */
|
} /* end loop age */
|
} /* end loop age */
|
fprintf(ficgp,"\nset out \"varpijgr%s%d%1d%1d-%1d%1d.png\";replot;",optionfilefiname, j1,k2,l2,k1,l1);
|
fprintf(ficgp,"\nset out \"varpijgr%s%d%1d%1d-%1d%1d.png\";replot;",optionfilefiname, j1,k1,l1,k2,l2);
|
first=1;
|
first=1;
|
} /*l12 */
|
} /*l12 */
|
} /* k12 */
|
} /* k12 */
|