version 1.11, 2001/05/17 16:07:14
|
version 1.12, 2002/02/20 16:57:00
|
Line 22
|
Line 22
|
delay between waves is not identical for each individual, or if some
|
delay between waves is not identical for each individual, or if some
|
individual missed an interview, the information is not rounded or lost, but
|
individual missed an interview, the information is not rounded or lost, but
|
taken into account using an interpolation or extrapolation.
|
taken into account using an interpolation or extrapolation.
|
hPijx is the probability to be
|
hPijx is the probability to be
|
observed in state i at age x+h conditional to the observed state i at age
|
observed in state i at age x+h conditional to the observed state i at age
|
x. The delay 'h' can be split into an exact number (nh*stepm) of
|
x. The delay 'h' can be split into an exact number (nh*stepm) of
|
unobserved intermediate states. This elementary transition (by month or
|
unobserved intermediate states. This elementary transition (by month or
|
Line 693 double **prevalim(double **prlim, int nl
|
Line 693 double **prevalim(double **prlim, int nl
|
}
|
}
|
}
|
}
|
|
|
/*************** transition probabilities **********/
|
/*************** transition probabilities ***************/
|
|
|
double **pmij(double **ps, double *cov, int ncovmodel, double *x, int nlstate )
|
double **pmij(double **ps, double *cov, int ncovmodel, double *x, int nlstate )
|
{
|
{
|
Line 716 double **pmij(double **ps, double *cov,
|
Line 716 double **pmij(double **ps, double *cov,
|
s2 += x[(i-1)*nlstate*ncovmodel+(j-2)*ncovmodel+nc+(i-1)*(ndeath-1)*ncovmodel]*cov[nc];
|
s2 += x[(i-1)*nlstate*ncovmodel+(j-2)*ncovmodel+nc+(i-1)*(ndeath-1)*ncovmodel]*cov[nc];
|
/*printf("Int j>i s1=%.17e, s2=%.17e %lx %lx\n",s1,s2,s1,s2);*/
|
/*printf("Int j>i s1=%.17e, s2=%.17e %lx %lx\n",s1,s2,s1,s2);*/
|
}
|
}
|
ps[i][j]=s2;
|
ps[i][j]=(s2);
|
}
|
}
|
}
|
}
|
|
/*ps[3][2]=1;*/
|
|
|
for(i=1; i<= nlstate; i++){
|
for(i=1; i<= nlstate; i++){
|
s1=0;
|
s1=0;
|
for(j=1; j<i; j++)
|
for(j=1; j<i; j++)
|
Line 740 double **pmij(double **ps, double *cov,
|
Line 742 double **pmij(double **ps, double *cov,
|
}
|
}
|
}
|
}
|
|
|
|
|
/* for(ii=1; ii<= nlstate+ndeath; ii++){
|
/* for(ii=1; ii<= nlstate+ndeath; ii++){
|
for(jj=1; jj<= nlstate+ndeath; jj++){
|
for(jj=1; jj<= nlstate+ndeath; jj++){
|
printf("%lf ",ps[ii][jj]);
|
printf("%lf ",ps[ii][jj]);
|
Line 804 double ***hpxij(double ***po, int nhstep
|
Line 807 double ***hpxij(double ***po, int nhstep
|
cov[1]=1.;
|
cov[1]=1.;
|
cov[2]=age+((h-1)*hstepm + (d-1))*stepm/YEARM;
|
cov[2]=age+((h-1)*hstepm + (d-1))*stepm/YEARM;
|
for (k=1; k<=cptcovn;k++) cov[2+k]=nbcode[Tvar[k]][codtab[ij][Tvar[k]]];
|
for (k=1; k<=cptcovn;k++) cov[2+k]=nbcode[Tvar[k]][codtab[ij][Tvar[k]]];
|
for (k=1; k<=cptcovage;k++)
|
for (k=1; k<=cptcovage;k++)
|
cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2];
|
cov[2+Tage[k]]=cov[2+Tage[k]]*cov[2];
|
for (k=1; k<=cptcovprod;k++)
|
for (k=1; k<=cptcovprod;k++)
|
cov[2+Tprod[k]]=nbcode[Tvard[k][1]][codtab[ij][Tvard[k][1]]]*nbcode[Tvard[k][2]][codtab[ij][Tvard[k][2]]];
|
cov[2+Tprod[k]]=nbcode[Tvard[k][1]][codtab[ij][Tvard[k][1]]]*nbcode[Tvard[k][2]][codtab[ij][Tvard[k][2]]];
|
|
|
|
|
Line 913 void hesscov(double **matcov, double p[]
|
Line 916 void hesscov(double **matcov, double p[]
|
void lubksb(double **a, int npar, int *indx, double b[]) ;
|
void lubksb(double **a, int npar, int *indx, double b[]) ;
|
void ludcmp(double **a, int npar, int *indx, double *d) ;
|
void ludcmp(double **a, int npar, int *indx, double *d) ;
|
|
|
|
|
hess=matrix(1,npar,1,npar);
|
hess=matrix(1,npar,1,npar);
|
|
|
printf("\nCalculation of the hessian matrix. Wait...\n");
|
printf("\nCalculation of the hessian matrix. Wait...\n");
|
Line 921 void hesscov(double **matcov, double p[]
|
Line 923 void hesscov(double **matcov, double p[]
|
printf("%d",i);fflush(stdout);
|
printf("%d",i);fflush(stdout);
|
hess[i][i]=hessii(p,ftolhess,i,delti);
|
hess[i][i]=hessii(p,ftolhess,i,delti);
|
/*printf(" %f ",p[i]);*/
|
/*printf(" %f ",p[i]);*/
|
|
/*printf(" %lf ",hess[i][i]);*/
|
}
|
}
|
|
|
for (i=1;i<=npar;i++) {
|
for (i=1;i<=npar;i++) {
|
for (j=1;j<=npar;j++) {
|
for (j=1;j<=npar;j++) {
|
if (j>i) {
|
if (j>i) {
|
printf(".%d%d",i,j);fflush(stdout);
|
printf(".%d%d",i,j);fflush(stdout);
|
hess[i][j]=hessij(p,delti,i,j);
|
hess[i][j]=hessij(p,delti,i,j);
|
hess[j][i]=hess[i][j];
|
hess[j][i]=hess[i][j];
|
|
/*printf(" %lf ",hess[i][j]);*/
|
}
|
}
|
}
|
}
|
}
|
}
|
Line 1032 double hessii( double x[], double delta,
|
Line 1036 double hessii( double x[], double delta,
|
}
|
}
|
}
|
}
|
delti[theta]=delts;
|
delti[theta]=delts;
|
return res;
|
return res;
|
|
|
}
|
}
|
|
|
Line 1311 void concatwav(int wav[], int **dh, int
|
Line 1315 void concatwav(int wav[], int **dh, int
|
if (j >= jmax) jmax=j;
|
if (j >= jmax) jmax=j;
|
if (j <= jmin) jmin=j;
|
if (j <= jmin) jmin=j;
|
sum=sum+j;
|
sum=sum+j;
|
if (j<0) printf("j=%d num=%d ",j,i);
|
/* if (j<10) printf("j=%d num=%d ",j,i); */
|
}
|
}
|
}
|
}
|
else{
|
else{
|
Line 1319 void concatwav(int wav[], int **dh, int
|
Line 1323 void concatwav(int wav[], int **dh, int
|
k=k+1;
|
k=k+1;
|
if (j >= jmax) jmax=j;
|
if (j >= jmax) jmax=j;
|
else if (j <= jmin)jmin=j;
|
else if (j <= jmin)jmin=j;
|
|
/* if (j<10) printf("j=%d jmin=%d num=%d ",j,jmin,i); */
|
sum=sum+j;
|
sum=sum+j;
|
}
|
}
|
jk= j/stepm;
|
jk= j/stepm;
|
Line 1335 void concatwav(int wav[], int **dh, int
|
Line 1340 void concatwav(int wav[], int **dh, int
|
}
|
}
|
jmean=sum/k;
|
jmean=sum/k;
|
printf("Delay (in months) between two waves Min=%d Max=%d Mean=%f\n\n ",jmin, jmax,jmean);
|
printf("Delay (in months) between two waves Min=%d Max=%d Mean=%f\n\n ",jmin, jmax,jmean);
|
}
|
}
|
/*********** Tricode ****************************/
|
/*********** Tricode ****************************/
|
void tricode(int *Tvar, int **nbcode, int imx)
|
void tricode(int *Tvar, int **nbcode, int imx)
|
{
|
{
|
Line 1373 void tricode(int *Tvar, int **nbcode, in
|
Line 1378 void tricode(int *Tvar, int **nbcode, in
|
|
|
for (k=0; k<19; k++) Ndum[k]=0;
|
for (k=0; k<19; k++) Ndum[k]=0;
|
|
|
for (i=1; i<=ncovmodel; i++) {
|
for (i=1; i<=ncovmodel-2; i++) {
|
ij=Tvar[i];
|
ij=Tvar[i];
|
Ndum[ij]++;
|
Ndum[ij]++;
|
}
|
}
|
Line 1448 void varevsij(char fileres[], double ***
|
Line 1453 void varevsij(char fileres[], double ***
|
double **dnewm,**doldm;
|
double **dnewm,**doldm;
|
int i, j, nhstepm, hstepm, h;
|
int i, j, nhstepm, hstepm, h;
|
int k, cptcode;
|
int k, cptcode;
|
double *xp;
|
double *xp;
|
double **gp, **gm;
|
double **gp, **gm;
|
double ***gradg, ***trgradg;
|
double ***gradg, ***trgradg;
|
double ***p3mat;
|
double ***p3mat;
|
Line 1663 int main()
|
Line 1668 int main()
|
int ju,jl, mi;
|
int ju,jl, mi;
|
int i1,j1, k1,k2,k3,jk,aa,bb, stepsize, ij;
|
int i1,j1, k1,k2,k3,jk,aa,bb, stepsize, ij;
|
int jnais,jdc,jint4,jint1,jint2,jint3,**outcome,**adl,*tab;
|
int jnais,jdc,jint4,jint1,jint2,jint3,**outcome,**adl,*tab;
|
|
|
int hstepm, nhstepm;
|
int hstepm, nhstepm;
|
double bage, fage, age, agelim, agebase;
|
double bage, fage, age, agelim, agebase;
|
double ftolpl=FTOL;
|
double ftolpl=FTOL;
|
Line 1772 split(pathtot, path,optionfile);
|
Line 1777 split(pathtot, path,optionfile);
|
fprintf(ficparo,"\n");
|
fprintf(ficparo,"\n");
|
}
|
}
|
|
|
npar= (nlstate+ndeath-1)*nlstate*ncovmodel;
|
npar= (nlstate+ndeath-1)*nlstate*ncovmodel;
|
|
|
p=param[1][1];
|
p=param[1][1];
|
|
|
/* Reads comments: lines beginning with '#' */
|
/* Reads comments: lines beginning with '#' */
|
Line 1862 split(pathtot, path,optionfile);
|
Line 1868 split(pathtot, path,optionfile);
|
tab=ivector(1,NCOVMAX);
|
tab=ivector(1,NCOVMAX);
|
ncodemax=ivector(1,8);
|
ncodemax=ivector(1,8);
|
|
|
i=1;
|
i=1;
|
while (fgets(line, MAXLINE, fic) != NULL) {
|
while (fgets(line, MAXLINE, fic) != NULL) {
|
if ((i >= firstobs) && (i <=lastobs)) {
|
if ((i >= firstobs) && (i <=lastobs)) {
|
|
|
Line 1884 split(pathtot, path,optionfile);
|
Line 1890 split(pathtot, path,optionfile);
|
cutv(stra, strb,line,' '); covar[j][i]=(double)(atoi(strb)); strcpy(line,stra);
|
cutv(stra, strb,line,' '); covar[j][i]=(double)(atoi(strb)); strcpy(line,stra);
|
}
|
}
|
num[i]=atol(stra);
|
num[i]=atol(stra);
|
|
|
/*printf("%d %.lf %.lf %.lf %.lf/%.lf %.lf/%.lf %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d\n",num[i],(covar[1][i]), (covar[2][i]), (weight[i]), (moisnais[i]), (annais[i]), (moisdc[i]), (andc[i]), (mint[1][i]), (anint[1][i]), (s[1][i]), (mint[2][i]), (anint[2][i]), (s[2][i]), (mint[3][i]), (anint[3][i]), (s[3][i]), (mint[4][i]), (anint[4][i]), (s[4][i]));*/
|
/*if((s[2][i]==2) && (s[3][i]==-1)&&(s[4][i]==9)){
|
|
printf("%d %.lf %.lf %.lf %.lf/%.lf %.lf/%.lf %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d\n",num[i],(covar[1][i]), (covar[2][i]),weight[i], (moisnais[i]), (annais[i]), (moisdc[i]), (andc[i]), (mint[1][i]), (anint[1][i]), (s[1][i]), (mint[2][i]), (anint[2][i]), (s[2][i]), (mint[3][i]), (anint[3][i]), (s[3][i]), (mint[4][i]), (anint[4][i]), (s[4][i])); ij=ij+1;}*/
|
|
|
i=i+1;
|
i=i+1;
|
}
|
}
|
}
|
}
|
|
/* printf("ii=%d", ij);
|
/*scanf("%d",i);*/
|
scanf("%d",i);*/
|
imx=i-1; /* Number of individuals */
|
imx=i-1; /* Number of individuals */
|
|
|
|
/* for (i=1; i<=imx; i++){
|
|
if ((s[1][i]==3) && (s[2][i]==2)) s[2][i]=3;
|
|
if ((s[2][i]==3) && (s[3][i]==2)) s[3][i]=3;
|
|
if ((s[3][i]==3) && (s[4][i]==2)) s[4][i]=3;
|
|
}
|
|
for (i=1; i<=imx; i++) printf("%d %.lf %.lf %.lf %.lf/%.lf %.lf/%.lf %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d %.lf/%.lf %d\n",num[i],(covar[1][i]), (covar[2][i]), (weight[i]), (moisnais[i]), (annais[i]), (moisdc[i]), (andc[i]), (mint[1][i]), (anint[1][i]), (s[1][i]), (mint[2][i]), (anint[2][i]), (s[2][i]), (mint[3][i]), (anint[3][i]), (s[3][i]), (mint[4][i]), (anint[4][i]), (s[4][i]));*/
|
|
|
/* Calculation of the number of parameter from char model*/
|
/* Calculation of the number of parameter from char model*/
|
Tvar=ivector(1,15);
|
Tvar=ivector(1,15);
|
Tprod=ivector(1,15);
|
Tprod=ivector(1,15);
|
Line 1975 split(pathtot, path,optionfile);
|
Line 1989 split(pathtot, path,optionfile);
|
}
|
}
|
/*-calculation of age at interview from date of interview and age at death -*/
|
/*-calculation of age at interview from date of interview and age at death -*/
|
agev=matrix(1,maxwav,1,imx);
|
agev=matrix(1,maxwav,1,imx);
|
|
|
|
for (i=1; i<=imx; i++)
|
|
for(m=2; (m<= maxwav); m++)
|
|
if ((mint[m][i]== 99) && (s[m][i] <= nlstate)){
|
|
anint[m][i]=9999;
|
|
s[m][i]=-1;
|
|
}
|
|
|
for (i=1; i<=imx; i++) {
|
for (i=1; i<=imx; i++) {
|
agedc[i]=(moisdc[i]/12.+andc[i])-(moisnais[i]/12.+annais[i]);
|
agedc[i]=(moisdc[i]/12.+andc[i])-(moisnais[i]/12.+annais[i]);
|
Line 2083 printf("Total number of individuals= %d,
|
Line 2104 printf("Total number of individuals= %d,
|
newms= matrix(1,nlstate+ndeath,1,nlstate+ndeath); /* creation */
|
newms= matrix(1,nlstate+ndeath,1,nlstate+ndeath); /* creation */
|
savms= matrix(1,nlstate+ndeath,1,nlstate+ndeath); /* creation */
|
savms= matrix(1,nlstate+ndeath,1,nlstate+ndeath); /* creation */
|
oldm=oldms; newm=newms; savm=savms; /* Keeps fixed addresses to free */
|
oldm=oldms; newm=newms; savm=savms; /* Keeps fixed addresses to free */
|
|
|
/* For Powell, parameters are in a vector p[] starting at p[1]
|
/* For Powell, parameters are in a vector p[] starting at p[1]
|
so we point p on param[1][1] so that p[1] maps on param[1][1][1] */
|
so we point p on param[1][1] so that p[1] maps on param[1][1][1] */
|
p=param[1][1]; /* *(*(*(param +1)+1)+0) */
|
p=param[1][1]; /* *(*(*(param +1)+1)+0) */
|
Line 2666 strcpy(fileresvpl,"vpl");
|
Line 2687 strcpy(fileresvpl,"vpl");
|
/*printf("Total time was %d uSec.\n", total_usecs);*/
|
/*printf("Total time was %d uSec.\n", total_usecs);*/
|
/*------ End -----------*/
|
/*------ End -----------*/
|
|
|
|
|
end:
|
end:
|
#ifdef windows
|
#ifdef windows
|
chdir(pathcd);
|
chdir(pathcd);
|