|
version 1.100, 2004/07/12 18:29:06
|
version 1.111, 2006/01/25 20:38:18
|
|
Line 1
|
Line 1
|
| /* $Id$ |
/* $Id$ |
| $State$ |
$State$ |
| $Log$ |
$Log$ |
| |
Revision 1.111 2006/01/25 20:38:18 brouard |
| |
(Module): Lots of cleaning and bugs added (Gompertz) |
| |
(Module): Comments can be added in data file. Missing date values |
| |
can be a simple dot '.'. |
| |
|
| |
Revision 1.110 2006/01/25 00:51:50 brouard |
| |
(Module): Lots of cleaning and bugs added (Gompertz) |
| |
|
| |
Revision 1.109 2006/01/24 19:37:15 brouard |
| |
(Module): Comments (lines starting with a #) are allowed in data. |
| |
|
| |
Revision 1.108 2006/01/19 18:05:42 lievre |
| |
Gnuplot problem appeared... |
| |
To be fixed |
| |
|
| |
Revision 1.107 2006/01/19 16:20:37 brouard |
| |
Test existence of gnuplot in imach path |
| |
|
| |
Revision 1.106 2006/01/19 13:24:36 brouard |
| |
Some cleaning and links added in html output |
| |
|
| |
Revision 1.105 2006/01/05 20:23:19 lievre |
| |
*** empty log message *** |
| |
|
| |
Revision 1.104 2005/09/30 16:11:43 lievre |
| |
(Module): sump fixed, loop imx fixed, and simplifications. |
| |
(Module): If the status is missing at the last wave but we know |
| |
that the person is alive, then we can code his/her status as -2 |
| |
(instead of missing=-1 in earlier versions) and his/her |
| |
contributions to the likelihood is 1 - Prob of dying from last |
| |
health status (= 1-p13= p11+p12 in the easiest case of somebody in |
| |
the healthy state at last known wave). Version is 0.98 |
| |
|
| |
Revision 1.103 2005/09/30 15:54:49 lievre |
| |
(Module): sump fixed, loop imx fixed, and simplifications. |
| |
|
| |
Revision 1.102 2004/09/15 17:31:30 brouard |
| |
Add the possibility to read data file including tab characters. |
| |
|
| |
Revision 1.101 2004/09/15 10:38:38 brouard |
| |
Fix on curr_time |
| |
|
| Revision 1.100 2004/07/12 18:29:06 brouard |
Revision 1.100 2004/07/12 18:29:06 brouard |
| Add version for Mac OS X. Just define UNIX in Makefile |
Add version for Mac OS X. Just define UNIX in Makefile |
| |
|
|
Line 200
|
Line 242
|
| #include <math.h> |
#include <math.h> |
| #include <stdio.h> |
#include <stdio.h> |
| #include <stdlib.h> |
#include <stdlib.h> |
| |
#include <string.h> |
| #include <unistd.h> |
#include <unistd.h> |
| |
|
| |
#include <limits.h> |
| |
#include <sys/types.h> |
| |
#include <sys/stat.h> |
| |
#include <errno.h> |
| |
extern int errno; |
| |
|
| /* #include <sys/time.h> */ |
/* #include <sys/time.h> */ |
| #include <time.h> |
#include <time.h> |
| #include "timeval.h" |
#include "timeval.h" |
|
Line 210
|
Line 259
|
| /* #define _(String) gettext (String) */ |
/* #define _(String) gettext (String) */ |
| |
|
| #define MAXLINE 256 |
#define MAXLINE 256 |
| |
|
| #define GNUPLOTPROGRAM "gnuplot" |
#define GNUPLOTPROGRAM "gnuplot" |
| /*#define GNUPLOTPROGRAM "..\\gp37mgw\\wgnuplot"*/ |
/*#define GNUPLOTPROGRAM "..\\gp37mgw\\wgnuplot"*/ |
| #define FILENAMELENGTH 132 |
#define FILENAMELENGTH 132 |
| /*#define DEBUG*/ |
|
| /*#define windows*/ |
|
| #define GLOCK_ERROR_NOPATH -1 /* empty path */ |
#define GLOCK_ERROR_NOPATH -1 /* empty path */ |
| #define GLOCK_ERROR_GETCWD -2 /* cannot get cwd */ |
#define GLOCK_ERROR_GETCWD -2 /* cannot get cwd */ |
| |
|
|
Line 232
|
Line 281
|
| #define AGEGOMP 10. /* Minimal age for Gompertz adjustment */ |
#define AGEGOMP 10. /* Minimal age for Gompertz adjustment */ |
| #ifdef UNIX |
#ifdef UNIX |
| #define DIRSEPARATOR '/' |
#define DIRSEPARATOR '/' |
| |
#define CHARSEPARATOR "/" |
| #define ODIRSEPARATOR '\\' |
#define ODIRSEPARATOR '\\' |
| #else |
#else |
| #define DIRSEPARATOR '\\' |
#define DIRSEPARATOR '\\' |
| |
#define CHARSEPARATOR "\\" |
| #define ODIRSEPARATOR '/' |
#define ODIRSEPARATOR '/' |
| #endif |
#endif |
| |
|
| /* $Id$ */ |
/* $Id$ */ |
| /* $State$ */ |
/* $State$ */ |
| |
|
| char version[]="Imach version 0.97b, May 2004, INED-EUROREVES "; |
char version[]="Imach version 0.98b, January 2006, INED-EUROREVES "; |
| char fullversion[]="$Revision$ $Date$"; |
char fullversion[]="$Revision$ $Date$"; |
| int erreur, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */ |
int erreur, nberr=0, nbwarn=0; /* Error number, number of errors number of warnings */ |
| int nvar; |
int nvar; |
|
Line 255 int popbased=0;
|
Line 306 int popbased=0;
|
| int *wav; /* Number of waves for this individuual 0 is possible */ |
int *wav; /* Number of waves for this individuual 0 is possible */ |
| int maxwav; /* Maxim number of waves */ |
int maxwav; /* Maxim number of waves */ |
| int jmin, jmax; /* min, max spacing between 2 waves */ |
int jmin, jmax; /* min, max spacing between 2 waves */ |
| |
int ijmin, ijmax; /* Individuals having jmin and jmax */ |
| int gipmx, gsw; /* Global variables on the number of contributions |
int gipmx, gsw; /* Global variables on the number of contributions |
| to the likelihood and the sum of weights (done by funcone)*/ |
to the likelihood and the sum of weights (done by funcone)*/ |
| int mle, weightopt; |
int mle, weightopt; |
|
Line 285 FILE *ficresvpl;
|
Line 337 FILE *ficresvpl;
|
| char fileresvpl[FILENAMELENGTH]; |
char fileresvpl[FILENAMELENGTH]; |
| char title[MAXLINE]; |
char title[MAXLINE]; |
| char optionfile[FILENAMELENGTH], datafile[FILENAMELENGTH], filerespl[FILENAMELENGTH]; |
char optionfile[FILENAMELENGTH], datafile[FILENAMELENGTH], filerespl[FILENAMELENGTH]; |
| char optionfilext[10], optionfilefiname[FILENAMELENGTH], plotcmd[FILENAMELENGTH]; |
char optionfilext[10], optionfilefiname[FILENAMELENGTH], plotcmd[FILENAMELENGTH], pplotcmd[FILENAMELENGTH]; |
| char tmpout[FILENAMELENGTH], tmpout2[FILENAMELENGTH]; |
char tmpout[FILENAMELENGTH], tmpout2[FILENAMELENGTH]; |
| char command[FILENAMELENGTH]; |
char command[FILENAMELENGTH]; |
| int outcmd=0; |
int outcmd=0; |
|
Line 307 long time_value;
|
Line 359 long time_value;
|
| extern long time(); |
extern long time(); |
| char strcurr[80], strfor[80]; |
char strcurr[80], strfor[80]; |
| |
|
| |
char *endptr; |
| |
long lval; |
| |
|
| #define NR_END 1 |
#define NR_END 1 |
| #define FREE_ARG char* |
#define FREE_ARG char* |
| #define FTOL 1.0e-10 |
#define FTOL 1.0e-10 |
|
Line 355 double *weight;
|
Line 410 double *weight;
|
| int **s; /* Status */ |
int **s; /* Status */ |
| double *agedc, **covar, idx; |
double *agedc, **covar, idx; |
| int **nbcode, *Tcode, *Tvar, **codtab, **Tvard, *Tprod, cptcovprod, *Tvaraff; |
int **nbcode, *Tcode, *Tvar, **codtab, **Tvard, *Tprod, cptcovprod, *Tvaraff; |
| |
double *lsurv, *lpop, *tpop; |
| |
|
| double ftol=FTOL; /* Tolerance for computing Max Likelihood */ |
double ftol=FTOL; /* Tolerance for computing Max Likelihood */ |
| double ftolhess; /* Tolerance for computing hessian */ |
double ftolhess; /* Tolerance for computing hessian */ |
|
Line 362 double ftolhess; /* Tolerance for comput
|
Line 418 double ftolhess; /* Tolerance for comput
|
| /**************** split *************************/ |
/**************** split *************************/ |
| static int split( char *path, char *dirc, char *name, char *ext, char *finame ) |
static int split( char *path, char *dirc, char *name, char *ext, char *finame ) |
| { |
{ |
| /* From a file name with full path (either Unix or Windows) we extract the directory (dirc) |
/* From a file name with (full) path (either Unix or Windows) we extract the directory (dirc) |
| the name of the file (name), its extension only (ext) and its first part of the name (finame) |
the name of the file (name), its extension only (ext) and its first part of the name (finame) |
| */ |
*/ |
| char *ss; /* pointer */ |
char *ss; /* pointer */ |
|
Line 371 static int split( char *path, char *dirc
|
Line 427 static int split( char *path, char *dirc
|
| l1 = strlen(path ); /* length of path */ |
l1 = strlen(path ); /* length of path */ |
| if ( l1 == 0 ) return( GLOCK_ERROR_NOPATH ); |
if ( l1 == 0 ) return( GLOCK_ERROR_NOPATH ); |
| ss= strrchr( path, DIRSEPARATOR ); /* find last / */ |
ss= strrchr( path, DIRSEPARATOR ); /* find last / */ |
| if ( ss == NULL ) { /* no directory, so use current */ |
if ( ss == NULL ) { /* no directory, so determine current directory */ |
| |
strcpy( name, path ); /* we got the fullname name because no directory */ |
| /*if(strrchr(path, ODIRSEPARATOR )==NULL) |
/*if(strrchr(path, ODIRSEPARATOR )==NULL) |
| printf("Warning you should use %s as a separator\n",DIRSEPARATOR);*/ |
printf("Warning you should use %s as a separator\n",DIRSEPARATOR);*/ |
| /* get current working directory */ |
/* get current working directory */ |
|
Line 379 static int split( char *path, char *dirc
|
Line 436 static int split( char *path, char *dirc
|
| if ( getcwd( dirc, FILENAME_MAX ) == NULL ) { |
if ( getcwd( dirc, FILENAME_MAX ) == NULL ) { |
| return( GLOCK_ERROR_GETCWD ); |
return( GLOCK_ERROR_GETCWD ); |
| } |
} |
| strcpy( name, path ); /* we've got it */ |
/* got dirc from getcwd*/ |
| |
printf(" DIRC = %s \n",dirc); |
| } else { /* strip direcotry from path */ |
} else { /* strip direcotry from path */ |
| ss++; /* after this, the filename */ |
ss++; /* after this, the filename */ |
| l2 = strlen( ss ); /* length of filename */ |
l2 = strlen( ss ); /* length of filename */ |
|
Line 387 static int split( char *path, char *dirc
|
Line 445 static int split( char *path, char *dirc
|
| strcpy( name, ss ); /* save file name */ |
strcpy( name, ss ); /* save file name */ |
| strncpy( dirc, path, l1 - l2 ); /* now the directory */ |
strncpy( dirc, path, l1 - l2 ); /* now the directory */ |
| dirc[l1-l2] = 0; /* add zero */ |
dirc[l1-l2] = 0; /* add zero */ |
| |
printf(" DIRC2 = %s \n",dirc); |
| } |
} |
| |
/* We add a separator at the end of dirc if not exists */ |
| l1 = strlen( dirc ); /* length of directory */ |
l1 = strlen( dirc ); /* length of directory */ |
| /*#ifdef windows |
if( dirc[l1-1] != DIRSEPARATOR ){ |
| if ( dirc[l1-1] != '\\' ) { dirc[l1] = '\\'; dirc[l1+1] = 0; } |
dirc[l1] = DIRSEPARATOR; |
| #else |
dirc[l1+1] = 0; |
| if ( dirc[l1-1] != '/' ) { dirc[l1] = '/'; dirc[l1+1] = 0; } |
printf(" DIRC3 = %s \n",dirc); |
| #endif |
} |
| */ |
|
| ss = strrchr( name, '.' ); /* find last / */ |
ss = strrchr( name, '.' ); /* find last / */ |
| if (ss >0){ |
if (ss >0){ |
| ss++; |
ss++; |
|
Line 404 static int split( char *path, char *dirc
|
Line 463 static int split( char *path, char *dirc
|
| strncpy( finame, name, l1-l2); |
strncpy( finame, name, l1-l2); |
| finame[l1-l2]= 0; |
finame[l1-l2]= 0; |
| } |
} |
| |
|
| return( 0 ); /* we're done */ |
return( 0 ); /* we're done */ |
| } |
} |
| |
|
|
Line 436 int nbocc(char *s, char occ)
|
Line 496 int nbocc(char *s, char occ)
|
| |
|
| void cutv(char *u,char *v, char*t, char occ) |
void cutv(char *u,char *v, char*t, char occ) |
| { |
{ |
| /* cuts string t into u and v where u is ended by char occ excluding it |
/* cuts string t into u and v where u ends before first occurence of char 'occ' |
| and v is after occ excluding it too : ex cutv(u,v,"abcdef2ghi2j",2) |
and v starts after first occurence of char 'occ' : ex cutv(u,v,"abcdef2ghi2j",'2') |
| gives u="abcedf" and v="ghi2j" */ |
gives u="abcedf" and v="ghi2j" */ |
| int i,lg,j,p=0; |
int i,lg,j,p=0; |
| i=0; |
i=0; |
|
Line 879 void powell(double p[], double **xi, int
|
Line 939 void powell(double p[], double **xi, int
|
| fprintf(ficrespow,"\n");fflush(ficrespow); |
fprintf(ficrespow,"\n");fflush(ficrespow); |
| if(*iter <=3){ |
if(*iter <=3){ |
| tm = *localtime(&curr_time.tv_sec); |
tm = *localtime(&curr_time.tv_sec); |
| strcpy(strcurr,asctime(&tmf)); |
strcpy(strcurr,asctime(&tm)); |
| /* asctime_r(&tm,strcurr); */ |
/* asctime_r(&tm,strcurr); */ |
| forecast_time=curr_time; |
forecast_time=curr_time; |
| itmp = strlen(strcurr); |
itmp = strlen(strcurr); |
| if(strcurr[itmp-1]=='\n') |
if(strcurr[itmp-1]=='\n') /* Windows outputs with a new line */ |
| strcurr[itmp-1]='\0'; |
strcurr[itmp-1]='\0'; |
| printf("\nConsidering the time needed for this last iteration #%d: %ld seconds,\n",*iter,curr_time.tv_sec-last_time.tv_sec); |
printf("\nConsidering the time needed for this last iteration #%d: %ld seconds,\n",*iter,curr_time.tv_sec-last_time.tv_sec); |
| fprintf(ficlog,"\nConsidering the time needed for this last iteration #%d: %ld seconds,\n",*iter,curr_time.tv_sec-last_time.tv_sec); |
fprintf(ficlog,"\nConsidering the time needed for this last iteration #%d: %ld seconds,\n",*iter,curr_time.tv_sec-last_time.tv_sec); |
|
Line 895 void powell(double p[], double **xi, int
|
Line 955 void powell(double p[], double **xi, int
|
| itmp = strlen(strfor); |
itmp = strlen(strfor); |
| if(strfor[itmp-1]=='\n') |
if(strfor[itmp-1]=='\n') |
| strfor[itmp-1]='\0'; |
strfor[itmp-1]='\0'; |
| printf(" - if your program needs %d iterations to converge, convergence will be \n reached in %s or\n on %s (current time is %s);\n",niterf, asc_diff_time(forecast_time.tv_sec-curr_time.tv_sec,tmpout),strfor,strcurr); |
printf(" - if your program needs %d iterations to converge, convergence will be \n reached in %s i.e.\n on %s (current time is %s);\n",niterf, asc_diff_time(forecast_time.tv_sec-curr_time.tv_sec,tmpout),strfor,strcurr); |
| fprintf(ficlog," - if your program needs %d iterations to converge, convergence will be \n reached in %s or\n on %s (current time is %s);\n",niterf, asc_diff_time(forecast_time.tv_sec-curr_time.tv_sec,tmpout),strfor,strcurr); |
fprintf(ficlog," - if your program needs %d iterations to converge, convergence will be \n reached in %s i.e.\n on %s (current time is %s);\n",niterf, asc_diff_time(forecast_time.tv_sec-curr_time.tv_sec,tmpout),strfor,strcurr); |
| } |
} |
| } |
} |
| for (i=1;i<=n;i++) { |
for (i=1;i<=n;i++) { |
|
Line 1235 double func( double *x)
|
Line 1295 double func( double *x)
|
| } /* end mult */ |
} /* end mult */ |
| |
|
| /*lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/ /* Original formula */ |
/*lli=log(out[s[mw[mi][i]][i]][s[mw[mi+1][i]][i]]);*/ /* Original formula */ |
| /* But now since version 0.9 we anticipate for bias and large stepm. |
/* But now since version 0.9 we anticipate for bias at large stepm. |
| * If stepm is larger than one month (smallest stepm) and if the exact delay |
* If stepm is larger than one month (smallest stepm) and if the exact delay |
| * (in months) between two waves is not a multiple of stepm, we rounded to |
* (in months) between two waves is not a multiple of stepm, we rounded to |
| * the nearest (and in case of equal distance, to the lowest) interval but now |
* the nearest (and in case of equal distance, to the lowest) interval but now |
| * we keep into memory the bias bh[mi][i] and also the previous matrix product |
* we keep into memory the bias bh[mi][i] and also the previous matrix product |
| * (i.e to dh[mi][i]-1) saved in 'savm'. The we inter(extra)polate the |
* (i.e to dh[mi][i]-1) saved in 'savm'. Then we inter(extra)polate the |
| * probability in order to take into account the bias as a fraction of the way |
* probability in order to take into account the bias as a fraction of the way |
| * from savm to out if bh is neagtive or even beyond if bh is positive. bh varies |
* from savm to out if bh is negative or even beyond if bh is positive. bh varies |
| * -stepm/2 to stepm/2 . |
* -stepm/2 to stepm/2 . |
| * For stepm=1 the results are the same as for previous versions of Imach. |
* For stepm=1 the results are the same as for previous versions of Imach. |
| * For stepm > 1 the results are less biased than in previous versions. |
* For stepm > 1 the results are less biased than in previous versions. |
|
Line 1250 double func( double *x)
|
Line 1310 double func( double *x)
|
| s1=s[mw[mi][i]][i]; |
s1=s[mw[mi][i]][i]; |
| s2=s[mw[mi+1][i]][i]; |
s2=s[mw[mi+1][i]][i]; |
| bbh=(double)bh[mi][i]/(double)stepm; |
bbh=(double)bh[mi][i]/(double)stepm; |
| /* bias is positive if real duration |
/* bias bh is positive if real duration |
| * is higher than the multiple of stepm and negative otherwise. |
* is higher than the multiple of stepm and negative otherwise. |
| */ |
*/ |
| /* lli= (savm[s1][s2]>1.e-8 ?(1.+bbh)*log(out[s1][s2])- bbh*log(savm[s1][s2]):log((1.+bbh)*out[s1][s2]));*/ |
/* lli= (savm[s1][s2]>1.e-8 ?(1.+bbh)*log(out[s1][s2])- bbh*log(savm[s1][s2]):log((1.+bbh)*out[s1][s2]));*/ |
| if( s2 > nlstate){ |
if( s2 > nlstate){ |
| /* i.e. if s2 is a death state and if the date of death is known then the contribution |
/* i.e. if s2 is a death state and if the date of death is known |
| to the likelihood is the probability to die between last step unit time and current |
then the contribution to the likelihood is the probability to |
| step unit time, which is also the differences between probability to die before dh |
die between last step unit time and current step unit time, |
| and probability to die before dh-stepm . |
which is also equal to probability to die before dh |
| |
minus probability to die before dh-stepm . |
| In version up to 0.92 likelihood was computed |
In version up to 0.92 likelihood was computed |
| as if date of death was unknown. Death was treated as any other |
as if date of death was unknown. Death was treated as any other |
| health state: the date of the interview describes the actual state |
health state: the date of the interview describes the actual state |
|
Line 1278 double func( double *x)
|
Line 1339 double func( double *x)
|
| lower mortality. |
lower mortality. |
| */ |
*/ |
| lli=log(out[s1][s2] - savm[s1][s2]); |
lli=log(out[s1][s2] - savm[s1][s2]); |
| }else{ |
|
| |
|
| |
} else if (s2==-2) { |
| |
for (j=1,survp=0. ; j<=nlstate; j++) |
| |
survp += out[s1][j]; |
| |
lli= survp; |
| |
} |
| |
|
| |
else if (s2==-4) { |
| |
for (j=3,survp=0. ; j<=nlstate; j++) |
| |
survp += out[s1][j]; |
| |
lli= survp; |
| |
} |
| |
|
| |
else if (s2==-5) { |
| |
for (j=1,survp=0. ; j<=2; j++) |
| |
survp += out[s1][j]; |
| |
lli= survp; |
| |
} |
| |
|
| |
|
| |
else{ |
| lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */ |
lli= log((1.+bbh)*out[s1][s2]- bbh*savm[s1][s2]); /* linear interpolation */ |
| /* lli= (savm[s1][s2]>(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*(savm[s1][s2])):log((1.+bbh)*out[s1][s2]));*/ /* linear interpolation */ |
/* lli= (savm[s1][s2]>(double)1.e-8 ?log((1.+bbh)*out[s1][s2]- bbh*(savm[s1][s2])):log((1.+bbh)*out[s1][s2]));*/ /* linear interpolation */ |
| } |
} |
|
Line 1842 void lubksb(double **a, int n, int *indx
|
Line 1924 void lubksb(double **a, int n, int *indx
|
| } |
} |
| |
|
| /************ Frequencies ********************/ |
/************ Frequencies ********************/ |
| void freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **agev, int nlstate, int imx, int *Tvaraff, int **nbcode, int *ncodemax,double **mint,double **anint) |
void freqsummary(char fileres[], int iagemin, int iagemax, int **s, double **agev, int nlstate, int imx, int *Tvaraff, int **nbcode, int *ncodemax,double **mint,double **anint, char strstart[]) |
| { /* Some frequencies */ |
{ /* Some frequencies */ |
| |
|
| int i, m, jk, k1,i1, j1, bool, z1,z2,j; |
int i, m, jk, k1,i1, j1, bool, z1,z2,j; |
|
Line 1862 void freqsummary(char fileres[], int ia
|
Line 1944 void freqsummary(char fileres[], int ia
|
| fprintf(ficlog,"Problem with prevalence resultfile: %s\n", fileresp); |
fprintf(ficlog,"Problem with prevalence resultfile: %s\n", fileresp); |
| exit(0); |
exit(0); |
| } |
} |
| freq= ma3x(-1,nlstate+ndeath,-1,nlstate+ndeath,iagemin,iagemax+3); |
freq= ma3x(-5,nlstate+ndeath,-5,nlstate+ndeath,iagemin,iagemax+3); |
| j1=0; |
j1=0; |
| |
|
| j=cptcoveff; |
j=cptcoveff; |
|
Line 1875 void freqsummary(char fileres[], int ia
|
Line 1957 void freqsummary(char fileres[], int ia
|
| j1++; |
j1++; |
| /*printf("cptcoveff=%d Tvaraff=%d", cptcoveff,Tvaraff[1]); |
/*printf("cptcoveff=%d Tvaraff=%d", cptcoveff,Tvaraff[1]); |
| scanf("%d", i);*/ |
scanf("%d", i);*/ |
| for (i=-1; i<=nlstate+ndeath; i++) |
for (i=-5; i<=nlstate+ndeath; i++) |
| for (jk=-1; jk<=nlstate+ndeath; jk++) |
for (jk=-5; jk<=nlstate+ndeath; jk++) |
| for(m=iagemin; m <= iagemax+3; m++) |
for(m=iagemin; m <= iagemax+3; m++) |
| freq[i][jk][m]=0; |
freq[i][jk][m]=0; |
| |
|
|
Line 1915 void freqsummary(char fileres[], int ia
|
Line 1997 void freqsummary(char fileres[], int ia
|
| } |
} |
| |
|
| /* fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/ |
/* fprintf(ficresp, "#Count between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf\n",jprev1, mprev1,anprev1,jprev2, mprev2,anprev2);*/ |
| |
fprintf(ficresp, "#Local time at start: %s", strstart); |
| if (cptcovn>0) { |
if (cptcovn>0) { |
| fprintf(ficresp, "\n#********** Variable "); |
fprintf(ficresp, "\n#********** Variable "); |
| for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresp, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtab[j1][z1]]); |
for (z1=1; z1<=cptcoveff; z1++) fprintf(ficresp, "V%d=%d ",Tvaraff[z1],nbcode[Tvaraff[z1]][codtab[j1][z1]]); |
|
Line 2001 void freqsummary(char fileres[], int ia
|
Line 2083 void freqsummary(char fileres[], int ia
|
| dateintmean=dateintsum/k2cpt; |
dateintmean=dateintsum/k2cpt; |
| |
|
| fclose(ficresp); |
fclose(ficresp); |
| free_ma3x(freq,-1,nlstate+ndeath,-1,nlstate+ndeath, iagemin, iagemax+3); |
free_ma3x(freq,-5,nlstate+ndeath,-5,nlstate+ndeath, iagemin, iagemax+3); |
| free_vector(pp,1,nlstate); |
free_vector(pp,1,nlstate); |
| free_matrix(prop,1,nlstate,iagemin, iagemax+3); |
free_matrix(prop,1,nlstate,iagemin, iagemax+3); |
| /* End of Freq */ |
/* End of Freq */ |
|
Line 2110 void concatwav(int wav[], int **dh, int
|
Line 2192 void concatwav(int wav[], int **dh, int
|
| mi=0; |
mi=0; |
| m=firstpass; |
m=firstpass; |
| while(s[m][i] <= nlstate){ |
while(s[m][i] <= nlstate){ |
| if(s[m][i]>=1) |
if(s[m][i]>=1 || s[m][i]==-2 || s[m][i]==-4 || s[m][i]==-5) |
| mw[++mi][i]=m; |
mw[++mi][i]=m; |
| if(m >=lastpass) |
if(m >=lastpass) |
| break; |
break; |
|
Line 2128 void concatwav(int wav[], int **dh, int
|
Line 2210 void concatwav(int wav[], int **dh, int
|
| if(mi==0){ |
if(mi==0){ |
| nbwarn++; |
nbwarn++; |
| if(first==0){ |
if(first==0){ |
| printf("Warning! None valid information for:%ld line=%d (skipped) and may be others, see log file\n",num[i],i); |
printf("Warning! No valid information for individual %ld line=%d (skipped) and may be others, see log file\n",num[i],i); |
| first=1; |
first=1; |
| } |
} |
| if(first==1){ |
if(first==1){ |
| fprintf(ficlog,"Warning! None valid information for:%ld line=%d (skipped)\n",num[i],i); |
fprintf(ficlog,"Warning! No valid information for individual %ld line=%d (skipped)\n",num[i],i); |
| } |
} |
| } /* end mi==0 */ |
} /* end mi==0 */ |
| } /* End individuals */ |
} /* End individuals */ |
|
Line 2150 void concatwav(int wav[], int **dh, int
|
Line 2232 void concatwav(int wav[], int **dh, int
|
| nberr++; |
nberr++; |
| printf("Error! Negative delay (%d to death) between waves %d and %d of individual %ld at line %d who is aged %.1f with statuses from %d to %d\n ",j,mw[mi][i],mw[mi+1][i],num[i], i,agev[mw[mi][i]][i],s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]); |
printf("Error! Negative delay (%d to death) between waves %d and %d of individual %ld at line %d who is aged %.1f with statuses from %d to %d\n ",j,mw[mi][i],mw[mi+1][i],num[i], i,agev[mw[mi][i]][i],s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]); |
| j=1; /* Temporary Dangerous patch */ |
j=1; /* Temporary Dangerous patch */ |
| printf(" We assumed that the date of interview was correct (and not the date of death) and postponed the death %d month(s) (one stepm) after the interview.\n You MUST fix the contradiction between dates.\n",stepm); |
printf(" We assumed that the date of interview was correct (and not the date of death) and postponed the death %d month(s) (one stepm) after the interview. You MUST fix the contradiction between dates.\n",stepm); |
| fprintf(ficlog,"Error! Negative delay (%d to death) between waves %d and %d of individual %ld at line %d who is aged %.1f with statuses from %d to %d\n ",j,mw[mi][i],mw[mi+1][i],num[i], i,agev[mw[mi][i]][i],s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]); |
fprintf(ficlog,"Error! Negative delay (%d to death) between waves %d and %d of individual %ld at line %d who is aged %.1f with statuses from %d to %d\n ",j,mw[mi][i],mw[mi+1][i],num[i], i,agev[mw[mi][i]][i],s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]); |
| fprintf(ficlog," We assumed that the date of interview was correct (and not the date of death) and postponed the death %d month(s) (one stepm) after the interview.\n You MUST fix the contradiction between dates.\n",stepm); |
fprintf(ficlog," We assumed that the date of interview was correct (and not the date of death) and postponed the death %d month(s) (one stepm) after the interview. You MUST fix the contradiction between dates.\n",stepm); |
| } |
} |
| k=k+1; |
k=k+1; |
| if (j >= jmax) jmax=j; |
if (j >= jmax){ |
| if (j <= jmin) jmin=j; |
jmax=j; |
| |
ijmax=i; |
| |
} |
| |
if (j <= jmin){ |
| |
jmin=j; |
| |
ijmin=i; |
| |
} |
| sum=sum+j; |
sum=sum+j; |
| /*if (j<0) printf("j=%d num=%d \n",j,i);*/ |
/*if (j<0) printf("j=%d num=%d \n",j,i);*/ |
| /* printf("%d %d %d %d\n", s[mw[mi][i]][i] ,s[mw[mi+1][i]][i],j,i);*/ |
/* printf("%d %d %d %d\n", s[mw[mi][i]][i] ,s[mw[mi+1][i]][i],j,i);*/ |
|
Line 2164 void concatwav(int wav[], int **dh, int
|
Line 2252 void concatwav(int wav[], int **dh, int
|
| } |
} |
| else{ |
else{ |
| j= rint( (agev[mw[mi+1][i]][i]*12 - agev[mw[mi][i]][i]*12)); |
j= rint( (agev[mw[mi+1][i]][i]*12 - agev[mw[mi][i]][i]*12)); |
| /* printf("%d %d %d %d\n", s[mw[mi][i]][i] ,s[mw[mi+1][i]][i],j,i);*/ |
/* if (j<0) printf("%d %lf %lf %d %d %d\n", i,agev[mw[mi+1][i]][i], agev[mw[mi][i]][i],j,s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]); */ |
| |
|
| k=k+1; |
k=k+1; |
| if (j >= jmax) jmax=j; |
if (j >= jmax) { |
| else if (j <= jmin)jmin=j; |
jmax=j; |
| |
ijmax=i; |
| |
} |
| |
else if (j <= jmin){ |
| |
jmin=j; |
| |
ijmin=i; |
| |
} |
| /* if (j<10) printf("j=%d jmin=%d num=%d ",j,jmin,i); */ |
/* if (j<10) printf("j=%d jmin=%d num=%d ",j,jmin,i); */ |
| /*printf("%d %lf %d %d %d\n", i,agev[mw[mi][i]][i],j,s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]);*/ |
/*printf("%d %lf %d %d %d\n", i,agev[mw[mi][i]][i],j,s[mw[mi][i]][i] ,s[mw[mi+1][i]][i]);*/ |
| if(j<0){ |
if(j<0){ |
|
Line 2210 void concatwav(int wav[], int **dh, int
|
Line 2305 void concatwav(int wav[], int **dh, int
|
| } /* end wave */ |
} /* end wave */ |
| } |
} |
| 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 (for indiviudal %ld) Max=%d (%ld) Mean=%f\n\n ",jmin, num[ijmin], jmax, num[ijmax], jmean); |
| fprintf(ficlog,"Delay (in months) between two waves Min=%d Max=%d Mean=%f\n\n ",jmin, jmax,jmean); |
fprintf(ficlog,"Delay (in months) between two waves Min=%d (for indiviudal %ld) Max=%d (%ld) Mean=%f\n\n ",jmin, ijmin, jmax, ijmax, jmean); |
| } |
} |
| |
|
| /*********** Tricode ****************************/ |
/*********** Tricode ****************************/ |
|
Line 2257 void tricode(int *Tvar, int **nbcode, in
|
Line 2352 void tricode(int *Tvar, int **nbcode, in
|
| for (k=0; k< maxncov; k++) Ndum[k]=0; |
for (k=0; k< maxncov; k++) Ndum[k]=0; |
| |
|
| for (i=1; i<=ncovmodel-2; i++) { |
for (i=1; i<=ncovmodel-2; i++) { |
| /* Listing of all covariables in staement model to see if some covariates appear twice. For example, V1 appears twice in V1+V1*V2.*/ |
/* Listing of all covariables in statement model to see if some covariates appear twice. For example, V1 appears twice in V1+V1*V2.*/ |
| ij=Tvar[i]; |
ij=Tvar[i]; |
| Ndum[ij]++; |
Ndum[ij]++; |
| } |
} |
|
Line 2275 void tricode(int *Tvar, int **nbcode, in
|
Line 2370 void tricode(int *Tvar, int **nbcode, in
|
| |
|
| /*********** Health Expectancies ****************/ |
/*********** Health Expectancies ****************/ |
| |
|
| void evsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int ij, int estepm,double delti[],double **matcov ) |
void evsij(char fileres[], double ***eij, double x[], int nlstate, int stepm, int bage, int fage, double **oldm, double **savm, int ij, int estepm,double delti[],double **matcov,char strstart[] ) |
| |
|
| { |
{ |
| /* Health expectancies */ |
/* Health expectancies */ |
|
Line 2293 void evsij(char fileres[], double ***eij
|
Line 2388 void evsij(char fileres[], double ***eij
|
| dnewm=matrix(1,nlstate*nlstate,1,npar); |
dnewm=matrix(1,nlstate*nlstate,1,npar); |
| doldm=matrix(1,nlstate*nlstate,1,nlstate*nlstate); |
doldm=matrix(1,nlstate*nlstate,1,nlstate*nlstate); |
| |
|
| |
fprintf(ficreseij,"# Local time at start: %s", strstart); |
| fprintf(ficreseij,"# Health expectancies\n"); |
fprintf(ficreseij,"# Health expectancies\n"); |
| fprintf(ficreseij,"# Age"); |
fprintf(ficreseij,"# Age"); |
| for(i=1; i<=nlstate;i++) |
for(i=1; i<=nlstate;i++) |
|
Line 2447 void evsij(char fileres[], double ***eij
|
Line 2543 void evsij(char fileres[], double ***eij
|
| } |
} |
| |
|
| /************ Variance ******************/ |
/************ Variance ******************/ |
| void varevsij(char optionfilefiname[], double ***vareij, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int ij, int estepm, int cptcov, int cptcod, int popbased, int mobilav) |
void varevsij(char optionfilefiname[], double ***vareij, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int ij, int estepm, int cptcov, int cptcod, int popbased, int mobilav, char strstart[]) |
| { |
{ |
| /* Variance of health expectancies */ |
/* Variance of health expectancies */ |
| /* double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double ** savm,double ftolpl);*/ |
/* double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double ** savm,double ftolpl);*/ |
|
Line 2498 void varevsij(char optionfilefiname[], d
|
Line 2594 void varevsij(char optionfilefiname[], d
|
| fprintf(ficlog,"Problem with resultfile: %s\n", fileresprobmorprev); |
fprintf(ficlog,"Problem with resultfile: %s\n", fileresprobmorprev); |
| } |
} |
| printf("Computing total mortality p.j=w1*p1j+w2*p2j+..: result on file '%s' \n",fileresprobmorprev); |
printf("Computing total mortality p.j=w1*p1j+w2*p2j+..: result on file '%s' \n",fileresprobmorprev); |
| |
|
| fprintf(ficlog,"Computing total mortality p.j=w1*p1j+w2*p2j+..: result on file '%s' \n",fileresprobmorprev); |
fprintf(ficlog,"Computing total mortality p.j=w1*p1j+w2*p2j+..: result on file '%s' \n",fileresprobmorprev); |
| |
fprintf(ficresprobmorprev, "#Local time at start: %s", strstart); |
| fprintf(ficresprobmorprev,"# probabilities of dying before estepm=%d months for people of exact age and weighted probabilities w1*p1j+w2*p2j+... stand dev in()\n",estepm); |
fprintf(ficresprobmorprev,"# probabilities of dying before estepm=%d months for people of exact age and weighted probabilities w1*p1j+w2*p2j+... stand dev in()\n",estepm); |
| fprintf(ficresprobmorprev,"# Age cov=%-d",ij); |
fprintf(ficresprobmorprev,"# Age cov=%-d",ij); |
| for(j=nlstate+1; j<=(nlstate+ndeath);j++){ |
for(j=nlstate+1; j<=(nlstate+ndeath);j++){ |
|
Line 2508 void varevsij(char optionfilefiname[], d
|
Line 2606 void varevsij(char optionfilefiname[], d
|
| } |
} |
| fprintf(ficresprobmorprev,"\n"); |
fprintf(ficresprobmorprev,"\n"); |
| fprintf(ficgp,"\n# Routine varevsij"); |
fprintf(ficgp,"\n# Routine varevsij"); |
| |
/* fprintf(fichtm, "#Local time at start: %s", strstart);*/ |
| fprintf(fichtm,"\n<li><h4> Computing probabilities of dying over estepm months as a weighted average (i.e global mortality independent of initial healh state)</h4></li>\n"); |
fprintf(fichtm,"\n<li><h4> Computing probabilities of dying over estepm months as a weighted average (i.e global mortality independent of initial healh state)</h4></li>\n"); |
| fprintf(fichtm,"\n<br>%s <br>\n",digitp); |
fprintf(fichtm,"\n<br>%s <br>\n",digitp); |
| /* } */ |
/* } */ |
| varppt = matrix(nlstate+1,nlstate+ndeath,nlstate+1,nlstate+ndeath); |
varppt = matrix(nlstate+1,nlstate+ndeath,nlstate+1,nlstate+ndeath); |
| |
fprintf(ficresvij, "#Local time at start: %s", strstart); |
| fprintf(ficresvij,"# Variance and covariance of health expectancies e.j \n# (weighted average of eij where weights are the stable prevalence in health states i\n"); |
fprintf(ficresvij,"# Variance and covariance of health expectancies e.j \n# (weighted average of eij where weights are the stable prevalence in health states i\n"); |
| fprintf(ficresvij,"# Age"); |
fprintf(ficresvij,"# Age"); |
| for(i=1; i<=nlstate;i++) |
for(i=1; i<=nlstate;i++) |
|
Line 2745 void varevsij(char optionfilefiname[], d
|
Line 2844 void varevsij(char optionfilefiname[], d
|
| } /* end varevsij */ |
} /* end varevsij */ |
| |
|
| /************ Variance of prevlim ******************/ |
/************ Variance of prevlim ******************/ |
| void varprevlim(char fileres[], double **varpl, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int ij) |
void varprevlim(char fileres[], double **varpl, double **matcov, double x[], double delti[], int nlstate, int stepm, double bage, double fage, double **oldm, double **savm, double **prlim, double ftolpl, int ij, char strstart[]) |
| { |
{ |
| /* Variance of prevalence limit */ |
/* Variance of prevalence limit */ |
| /* double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double **savm,double ftolpl);*/ |
/* double **prevalim(double **prlim, int nlstate, double *xp, double age, double **oldm, double **savm,double ftolpl);*/ |
|
Line 2758 void varprevlim(char fileres[], double *
|
Line 2857 void varprevlim(char fileres[], double *
|
| double **gradg, **trgradg; |
double **gradg, **trgradg; |
| double age,agelim; |
double age,agelim; |
| int theta; |
int theta; |
| |
fprintf(ficresvpl, "#Local time at start: %s", strstart); |
| fprintf(ficresvpl,"# Standard deviation of stable prevalences \n"); |
fprintf(ficresvpl,"# Standard deviation of stable prevalences \n"); |
| fprintf(ficresvpl,"# Age"); |
fprintf(ficresvpl,"# Age"); |
| for(i=1; i<=nlstate;i++) |
for(i=1; i<=nlstate;i++) |
|
Line 2828 void varprevlim(char fileres[], double *
|
Line 2927 void varprevlim(char fileres[], double *
|
| } |
} |
| |
|
| /************ 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) |
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, i1, k1, l1, t, tj; |
int i, j=0, i1, k1, l1, t, tj; |
| int k2, l2, j1, z1; |
int k2, l2, j1, z1; |
|
Line 2873 void varprob(char optionfilefiname[], do
|
Line 2972 void varprob(char optionfilefiname[], do
|
| 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); |
| |
fprintf(ficresprob, "#Local time at start: %s", strstart); |
| 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"); |
| |
fprintf(ficresprobcov, "#Local time at start: %s", strstart); |
| fprintf(ficresprobcov,"#One-step probabilities and covariance matrix\n"); |
fprintf(ficresprobcov,"#One-step probabilities and covariance matrix\n"); |
| fprintf(ficresprobcov,"# Age"); |
fprintf(ficresprobcov,"# Age"); |
| |
fprintf(ficresprobcor, "#Local time at start: %s", strstart); |
| fprintf(ficresprobcor,"#One-step probabilities and correlation matrix\n"); |
fprintf(ficresprobcor,"#One-step probabilities and correlation matrix\n"); |
| fprintf(ficresprobcov,"# Age"); |
fprintf(ficresprobcov,"# Age"); |
| |
|
|
Line 3148 void printinghtml(char fileres[], char t
|
Line 3249 void printinghtml(char fileres[], char t
|
| double jprev2, double mprev2,double anprev2){ |
double jprev2, double mprev2,double anprev2){ |
| int jj1, k1, i1, cpt; |
int jj1, k1, i1, cpt; |
| |
|
| fprintf(fichtm,"<ul><li><h4>Result files (first order: no variance)</h4>\n \ |
fprintf(fichtm,"<ul><li><a href='#firstorder'>Result files (first order: no variance)</a>\n \ |
| |
<li><a href='#secondorder'>Result files (second order (variance)</a>\n \ |
| |
</ul>"); |
| |
fprintf(fichtm,"<ul><li><h4><a name='firstorder'>Result files (first order: no variance)</a></h4>\n \ |
| - Observed prevalence in each state (during the period defined between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf): <a href=\"%s\">%s</a> <br>\n ", |
- Observed prevalence in each state (during the period defined between %.lf/%.lf/%.lf and %.lf/%.lf/%.lf): <a href=\"%s\">%s</a> <br>\n ", |
| jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,subdirf2(fileres,"p"),subdirf2(fileres,"p")); |
jprev1, mprev1,anprev1,jprev2, mprev2,anprev2,subdirf2(fileres,"p"),subdirf2(fileres,"p")); |
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
|
Line 3199 fprintf(fichtm," \n<ul><li><b>Graphs</b>
|
Line 3303 fprintf(fichtm," \n<ul><li><b>Graphs</b>
|
| |
|
| |
|
| fprintf(fichtm,"\ |
fprintf(fichtm,"\ |
| \n<br><li><h4> Result files (second order: variances)</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>\n", rfileres,rfileres); |
- Parameter file with estimated parameters and covariance matrix: <a href=\"%s\">%s</a> <br>\n", rfileres,rfileres); |
| |
|
| fprintf(fichtm," - Variance of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
fprintf(fichtm," - Variance of one-step probabilities: <a href=\"%s\">%s</a> <br>\n", |
|
Line 3913 double gompertz(double x[])
|
Line 4017 double gompertz(double x[])
|
| { |
{ |
| double A,B,L=0.0,sump=0.,num=0.; |
double A,B,L=0.0,sump=0.,num=0.; |
| int i,n=0; /* n is the size of the sample */ |
int i,n=0; /* n is the size of the sample */ |
| |
|
| for (i=0;i<=imx-1 ; i++) { |
for (i=0;i<=imx-1 ; i++) { |
| sump=sump+weight[i]; |
sump=sump+weight[i]; |
| sump=sump+1; |
/* sump=sump+1;*/ |
| num=num+1; |
num=num+1; |
| } |
} |
| |
|
| |
|
| /* for (i=1; i<=imx; i++) |
/* for (i=0; i<=imx; i++) |
| if (wav[i]>0) printf("i=%d ageex=%lf agecens=%lf agedc=%lf cens=%d %d\n" ,i,ageexmed[i],agecens[i],agedc[i],cens[i],wav[i]);*/ |
if (wav[i]>0) printf("i=%d ageex=%lf agecens=%lf agedc=%lf cens=%d %d\n" ,i,ageexmed[i],agecens[i],agedc[i],cens[i],wav[i]);*/ |
| |
|
| for (i=0;i<=imx-1 ; i++) |
for (i=1;i<=imx ; i++) |
| { |
{ |
| if (cens[i]==1 & wav[i]>1) |
if (cens[i] == 1 && wav[i]>1) |
| A=-x[1]/(x[2])* |
A=-x[1]/(x[2])*(exp(x[2]*(agecens[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp))); |
| (exp(x[2]/YEARM*(agecens[i]*12-agegomp*12))-exp(x[2]/YEARM*(ageexmed[i]*12-agegomp*12))); |
|
| |
|
| if (cens[i]==0 & wav[i]>1) |
if (cens[i] == 0 && wav[i]>1) |
| A=-x[1]/(x[2])* |
A=-x[1]/(x[2])*(exp(x[2]*(agedc[i]-agegomp))-exp(x[2]*(ageexmed[i]-agegomp))) |
| (exp(x[2]/YEARM*(agedc[i]*12-agegomp*12))-exp(x[2]/YEARM*(ageexmed[i]*12-agegomp*12))) |
+log(x[1]/YEARM)+x[2]*(agedc[i]-agegomp)+log(YEARM); |
| +log(x[1]/YEARM)+x[2]/YEARM*(agedc[i]*12-agegomp*12)+log(YEARM); |
|
| |
|
| if (wav[i]>1 & agecens[i]>15) { |
/*if (wav[i] > 1 && agecens[i] > 15) {*/ /* ??? */ |
| |
if (wav[i] > 1 ) { /* ??? */ |
| L=L+A*weight[i]; |
L=L+A*weight[i]; |
| /* printf("\ni=%d A=%f L=%lf x[1]=%lf x[2]=%lf ageex=%lf agecens=%lf cens=%d agedc=%lf weight=%lf\n",i,A,L,x[1],x[2],ageexmed[i]*12,agecens[i]*12,cens[i],agedc[i]*12,weight[i]);*/ |
/* printf("\ni=%d A=%f L=%lf x[1]=%lf x[2]=%lf ageex=%lf agecens=%lf cens=%d agedc=%lf weight=%lf\n",i,A,L,x[1],x[2],ageexmed[i]*12,agecens[i]*12,cens[i],agedc[i]*12,weight[i]);*/ |
| } |
} |
|
Line 3948 double gompertz(double x[])
|
Line 4052 double gompertz(double x[])
|
| /******************* Printing html file ***********/ |
/******************* Printing html file ***********/ |
| void printinghtmlmort(char fileres[], char title[], char datafile[], int firstpass, \ |
void printinghtmlmort(char fileres[], char title[], char datafile[], int firstpass, \ |
| int lastpass, int stepm, int weightopt, char model[],\ |
int lastpass, int stepm, int weightopt, char model[],\ |
| int imx, double p[],double **matcov){ |
int imx, double p[],double **matcov,double agemortsup){ |
| int i; |
int i,k; |
| |
|
| fprintf(fichtm,"<ul><li><h4>Result files </h4>\n Force of mortality. Parameters of the Gompertz fit (with confidence interval in brackets):<br>"); |
fprintf(fichtm,"<ul><li><h4>Result files </h4>\n Force of mortality. Parameters of the Gompertz fit (with confidence interval in brackets):<br>"); |
| fprintf(fichtm," mu(age) =%lf*exp(%lf*(age-%d)) per year<br><br>",p[1],p[2],agegomp); |
fprintf(fichtm," mu(age) =%lf*exp(%lf*(age-%d)) per year<br><br>",p[1],p[2],agegomp); |
|
Line 3957 void printinghtmlmort(char fileres[], ch
|
Line 4061 void printinghtmlmort(char fileres[], ch
|
| fprintf(fichtm," p[%d] = %lf [%f ; %f]<br>\n",i,p[i],p[i]-2*sqrt(matcov[i][i]),p[i]+2*sqrt(matcov[i][i])); |
fprintf(fichtm," p[%d] = %lf [%f ; %f]<br>\n",i,p[i],p[i]-2*sqrt(matcov[i][i]),p[i]+2*sqrt(matcov[i][i])); |
| fprintf(fichtm,"<br><br><img src=\"graphmort.png\">"); |
fprintf(fichtm,"<br><br><img src=\"graphmort.png\">"); |
| fprintf(fichtm,"</ul>"); |
fprintf(fichtm,"</ul>"); |
| |
|
| |
fprintf(fichtm,"<ul><li><h4>Life table</h4>\n <br>"); |
| |
|
| |
fprintf(fichtm,"\nAge l<inf>x</inf> q<inf>x</inf> d(x,x+1) L<inf>x</inf> T<inf>x</inf> e<infx</inf><br>"); |
| |
|
| |
for (k=agegomp;k<(agemortsup-2);k++) |
| |
fprintf(fichtm,"%d %.0lf %lf %.0lf %.0lf %.0lf %lf<br>\n",k,lsurv[k],p[1]*exp(p[2]*(k-agegomp)),(p[1]*exp(p[2]*(k-agegomp)))*lsurv[k],lpop[k],tpop[k],tpop[k]/lsurv[k]); |
| |
|
| |
|
| fflush(fichtm); |
fflush(fichtm); |
| } |
} |
| |
|
|
Line 3994 int main(int argc, char *argv[])
|
Line 4107 int main(int argc, char *argv[])
|
| { |
{ |
| 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,j, k, n=MAXN,iter,m,size=100,cptcode, cptcod; |
int i,j, k, n=MAXN,iter,m,size=100,cptcode, cptcod; |
| |
int linei, month, year,iout; |
| int jj, ll, li, lj, lk, imk; |
int jj, ll, li, lj, lk, imk; |
| int numlinepar=0; /* Current linenumber of parameter file */ |
int numlinepar=0; /* Current linenumber of parameter file */ |
| int itimes; |
int itimes; |
| int NDIM=2; |
int NDIM=2; |
| |
|
| char ca[32], cb[32], cc[32]; |
char ca[32], cb[32], cc[32]; |
| |
char dummy[]=" "; |
| /* FILE *fichtm; *//* Html File */ |
/* FILE *fichtm; *//* Html File */ |
| /* FILE *ficgp;*/ /*Gnuplot File */ |
/* FILE *ficgp;*/ /*Gnuplot File */ |
| |
struct stat info; |
| double agedeb, agefin,hf; |
double agedeb, agefin,hf; |
| double ageminpar=1.e20,agemin=1.e20, agemaxpar=-1.e20, agemax=-1.e20; |
double ageminpar=1.e20,agemin=1.e20, agemaxpar=-1.e20, agemax=-1.e20; |
| |
|
|
Line 4024 int main(int argc, char *argv[])
|
Line 4140 int main(int argc, char *argv[])
|
| int mobilavproj=0 , prevfcast=0 ; /* moving average of prev, If prevfcast=1 prevalence projection */ |
int mobilavproj=0 , prevfcast=0 ; /* moving average of prev, If prevfcast=1 prevalence projection */ |
| int mobilav=0,popforecast=0; |
int mobilav=0,popforecast=0; |
| int hstepm, nhstepm; |
int hstepm, nhstepm; |
| |
int agemortsup; |
| |
float sumlpop=0.; |
| double jprev1=1, mprev1=1,anprev1=2000,jprev2=1, mprev2=1,anprev2=2000; |
double jprev1=1, mprev1=1,anprev1=2000,jprev2=1, mprev2=1,anprev2=2000; |
| double jpyram=1, mpyram=1,anpyram=2000,jpyram1=1, mpyram1=1,anpyram1=2000; |
double jpyram=1, mpyram=1,anpyram=2000,jpyram1=1, mpyram1=1,anpyram1=2000; |
| |
|
|
Line 4100 int main(int argc, char *argv[])
|
Line 4218 int main(int argc, char *argv[])
|
| printf("pathtot=%s, path=%s, optionfile=%s\n",pathtot,path,optionfile);*/ |
printf("pathtot=%s, path=%s, optionfile=%s\n",pathtot,path,optionfile);*/ |
| /* cutv(path,optionfile,pathtot,'\\');*/ |
/* cutv(path,optionfile,pathtot,'\\');*/ |
| |
|
| |
/* Split argv[0], imach program to get pathimach */ |
| |
printf("\nargv[0]=%s argv[1]=%s, \n",argv[0],argv[1]); |
| split(argv[0],pathimach,optionfile,optionfilext,optionfilefiname); |
split(argv[0],pathimach,optionfile,optionfilext,optionfilefiname); |
| |
printf("\nargv[0]=%s pathimach=%s, \noptionfile=%s \noptionfilext=%s \noptionfilefiname=%s\n",argv[0],pathimach,optionfile,optionfilext,optionfilefiname); |
| /* strcpy(pathimach,argv[0]); */ |
/* strcpy(pathimach,argv[0]); */ |
| |
/* Split argv[1]=pathtot, parameter file name to get path, optionfile, extension and name */ |
| split(pathtot,path,optionfile,optionfilext,optionfilefiname); |
split(pathtot,path,optionfile,optionfilext,optionfilefiname); |
| printf("pathimach=%s, pathtot=%s,\npath=%s,\noptionfile=%s \noptionfilext=%s \noptionfilefiname=%s\n",pathimach,pathtot,path,optionfile,optionfilext,optionfilefiname); |
printf("\npathtot=%s,\npath=%s,\noptionfile=%s \noptionfilext=%s \noptionfilefiname=%s\n",pathtot,path,optionfile,optionfilext,optionfilefiname); |
| chdir(path); |
chdir(path); |
| strcpy(command,"mkdir "); |
strcpy(command,"mkdir "); |
| strcat(command,optionfilefiname); |
strcat(command,optionfilefiname); |
|
Line 4270 int main(int argc, char *argv[])
|
Line 4392 int main(int argc, char *argv[])
|
| } |
} |
| fflush(ficlog); |
fflush(ficlog); |
| |
|
| |
|
| p=param[1][1]; |
p=param[1][1]; |
| |
|
| /* Reads comments: lines beginning with '#' */ |
/* Reads comments: lines beginning with '#' */ |
|
Line 4394 int main(int argc, char *argv[])
|
Line 4515 int main(int argc, char *argv[])
|
| ncodemax=ivector(1,8); |
ncodemax=ivector(1,8); |
| |
|
| i=1; |
i=1; |
| while (fgets(line, MAXLINE, fic) != NULL) { |
linei=0; |
| if ((i >= firstobs) && (i <=lastobs)) { |
while ((fgets(line, MAXLINE, fic) != NULL) &&((i >= firstobs) && (i <=lastobs))) { |
| |
linei=linei+1; |
| for (j=maxwav;j>=1;j--){ |
for(j=strlen(line); j>=0;j--){ /* Untabifies line */ |
| cutv(stra, strb,line,' '); s[j][i]=atoi(strb); |
if(line[j] == '\t') |
| strcpy(line,stra); |
line[j] = ' '; |
| cutv(stra, strb,line,'/'); anint[j][i]=(double)(atoi(strb)); strcpy(line,stra); |
} |
| cutv(stra, strb,line,' '); mint[j][i]=(double)(atoi(strb)); strcpy(line,stra); |
for(j=strlen(line)-1; (line[j]==' ')||(line[j]==10)||(line[j]==13);j--){ |
| |
; |
| |
}; |
| |
line[j+1]=0; /* Trims blanks at end of line */ |
| |
if(line[0]=='#'){ |
| |
fprintf(ficlog,"Comment line\n%s\n",line); |
| |
printf("Comment line\n%s\n",line); |
| |
continue; |
| |
} |
| |
|
| |
for (j=maxwav;j>=1;j--){ |
| |
cutv(stra, strb,line,' '); |
| |
errno=0; |
| |
lval=strtol(strb,&endptr,10); |
| |
/* if (errno == ERANGE && (lval == LONG_MAX || lval == LONG_MIN))*/ |
| |
if( strb[0]=='\0' || (*endptr != '\0')){ |
| |
printf("Error reading data around '%d' at line number %d %s for individual %d, '%s'\nShould be a status of wave %d. Setting maxwav=%d might be wrong. Exiting.\n", strb, linei,i,line,j,maxwav); |
| |
exit(1); |
| } |
} |
| |
s[j][i]=lval; |
| cutv(stra, strb,line,'/'); andc[i]=(double)(atoi(strb)); strcpy(line,stra); |
|
| cutv(stra, strb,line,' '); moisdc[i]=(double)(atoi(strb)); strcpy(line,stra); |
strcpy(line,stra); |
| |
cutv(stra, strb,line,' '); |
| cutv(stra, strb,line,'/'); annais[i]=(double)(atoi(strb)); strcpy(line,stra); |
if(iout=sscanf(strb,"%d/%d",&month, &year) != 0){ |
| cutv(stra, strb,line,' '); moisnais[i]=(double)(atoi(strb)); strcpy(line,stra); |
} |
| |
else if(iout=sscanf(strb,"%s.") != 0){ |
| cutv(stra, strb,line,' '); weight[i]=(double)(atoi(strb)); strcpy(line,stra); |
month=99; |
| for (j=ncovcol;j>=1;j--){ |
year=9999; |
| cutv(stra, strb,line,' '); covar[j][i]=(double)(atoi(strb)); strcpy(line,stra); |
}else{ |
| } |
printf("Error reading data around '%s' at line number %ld %s for individual %d, '%s'\nShould be a date of interview (mm/yyyy or .) at wave %d. Exiting.\n",strb, linei,i, line,j); |
| lstra=strlen(stra); |
exit(1); |
| if(lstra > 9){ /* More than 2**32 or max of what printf can write with %ld */ |
} |
| stratrunc = &(stra[lstra-9]); |
anint[j][i]= (double) year; |
| num[i]=atol(stratrunc); |
mint[j][i]= (double)month; |
| |
strcpy(line,stra); |
| |
} /* ENd Waves */ |
| |
|
| |
cutv(stra, strb,line,' '); |
| |
if(iout=sscanf(strb,"%d/%d",&month, &year) != 0){ |
| |
} |
| |
else if(iout=sscanf(strb,"%s.",dummy) != 0){ |
| |
month=99; |
| |
year=9999; |
| |
}else{ |
| |
printf("Error reading data around '%s' at line number %ld %s for individual %d, '%s'\nShould be a date of death (mm/yyyy or .). Exiting.\n",strb, linei,i,line); |
| |
exit(1); |
| |
} |
| |
andc[i]=(double) year; |
| |
moisdc[i]=(double) month; |
| |
strcpy(line,stra); |
| |
|
| |
cutv(stra, strb,line,' '); |
| |
if(iout=sscanf(strb,"%d/%d",&month, &year) != 0){ |
| |
} |
| |
else if(iout=sscanf(strb,"%s.") != 0){ |
| |
month=99; |
| |
year=9999; |
| |
}else{ |
| |
printf("Error reading data around '%s' at line number %ld %s for individual %d, '%s'\nShould be a date of birth (mm/yyyy or .). Exiting.\n",strb, linei,i,line,j); |
| |
exit(1); |
| |
} |
| |
annais[i]=(double)(year); |
| |
moisnais[i]=(double)(month); |
| |
strcpy(line,stra); |
| |
|
| |
cutv(stra, strb,line,' '); |
| |
errno=0; |
| |
lval=strtol(strb,&endptr,10); |
| |
if( strb[0]=='\0' || (*endptr != '\0')){ |
| |
printf("Error reading data around '%d' at line number %ld %s for individual %d\nShould be a weight. Exiting.\n",lval, i,line,linei); |
| |
exit(1); |
| |
} |
| |
weight[i]=(double)(lval); |
| |
strcpy(line,stra); |
| |
|
| |
for (j=ncovcol;j>=1;j--){ |
| |
cutv(stra, strb,line,' '); |
| |
errno=0; |
| |
lval=strtol(strb,&endptr,10); |
| |
if( strb[0]=='\0' || (*endptr != '\0')){ |
| |
printf("Error reading data around '%d' at line number %ld %s for individual %d, '%s'\nShould be a covar (meaning 0 for the reference or 1). Exiting.\n",lval, linei,i, line); |
| |
exit(1); |
| |
} |
| |
if(lval <-1 || lval >1){ |
| |
printf("Error reading data around '%d' at line number %ld %s for individual %d, '%s'\nShould be a value of the %d covar (meaning 0 for the reference or 1. IMaCh does not build design variables, do it your self). Exiting.\n",lval,linei, i,line,j); |
| |
exit(1); |
| } |
} |
| else |
covar[j][i]=(double)(lval); |
| num[i]=atol(stra); |
strcpy(line,stra); |
| |
} |
| /*if((s[2][i]==2) && (s[3][i]==-1)&&(s[4][i]==9)){ |
lstra=strlen(stra); |
| printf("%ld %.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;}*/ |
|
| |
if(lstra > 9){ /* More than 2**32 or max of what printf can write with %ld */ |
| i=i+1; |
stratrunc = &(stra[lstra-9]); |
| } |
num[i]=atol(stratrunc); |
| } |
} |
| |
else |
| |
num[i]=atol(stra); |
| |
/*if((s[2][i]==2) && (s[3][i]==-1)&&(s[4][i]==9)){ |
| |
printf("%ld %.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; |
| |
} /* End loop reading data */ |
| /* printf("ii=%d", ij); |
/* printf("ii=%d", ij); |
| scanf("%d",i);*/ |
scanf("%d",i);*/ |
| imx=i-1; /* Number of individuals */ |
imx=i-1; /* Number of individuals */ |
|
Line 4441 int main(int argc, char *argv[])
|
Line 4638 int main(int argc, char *argv[])
|
| if (s[4][i]==9) s[4][i]=-1; |
if (s[4][i]==9) s[4][i]=-1; |
| printf("%ld %.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]));}*/ |
printf("%ld %.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]));}*/ |
| |
|
| for (i=1; i<=imx; i++) |
/* for (i=1; i<=imx; i++) */ |
| |
|
| /*if ((s[3][i]==3) || (s[4][i]==3)) weight[i]=0.08; |
/*if ((s[3][i]==3) || (s[4][i]==3)) weight[i]=0.08; |
| else weight[i]=1;*/ |
else weight[i]=1;*/ |
| |
|
| /* Calculation of the number of parameter from char model*/ |
/* Calculation of the number of parameters from char model */ |
| Tvar=ivector(1,15); /* stores the number n of the covariates in Vm+Vn at 1 and m at 2 */ |
Tvar=ivector(1,15); /* stores the number n of the covariates in Vm+Vn at 1 and m at 2 */ |
| Tprod=ivector(1,15); |
Tprod=ivector(1,15); |
| Tvaraff=ivector(1,15); |
Tvaraff=ivector(1,15); |
|
Line 4559 int main(int argc, char *argv[])
|
Line 4756 int main(int argc, char *argv[])
|
| 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]); |
| for(m=firstpass; (m<= lastpass); m++){ |
for(m=firstpass; (m<= lastpass); m++){ |
| if(s[m][i] >0){ |
if(s[m][i] >0 || s[m][i]==-2 || s[m][i]==-4 || s[m][i]==-5){ |
| if (s[m][i] >= nlstate+1) { |
if (s[m][i] >= nlstate+1) { |
| if(agedc[i]>0) |
if(agedc[i]>0) |
| if((int)moisdc[i]!=99 && (int)andc[i]!=9999) |
if((int)moisdc[i]!=99 && (int)andc[i]!=9999) |
|
Line 4575 int main(int argc, char *argv[])
|
Line 4772 int main(int argc, char *argv[])
|
| } |
} |
| } |
} |
| else if(s[m][i] !=9){ /* Standard case, age in fractional |
else if(s[m][i] !=9){ /* Standard case, age in fractional |
| years but with the precision of a |
years but with the precision of a month */ |
| month */ |
|
| agev[m][i]=(mint[m][i]/12.+1./24.+anint[m][i])-(moisnais[i]/12.+1./24.+annais[i]); |
agev[m][i]=(mint[m][i]/12.+1./24.+anint[m][i])-(moisnais[i]/12.+1./24.+annais[i]); |
| if((int)mint[m][i]==99 || (int)anint[m][i]==9999) |
if((int)mint[m][i]==99 || (int)anint[m][i]==9999) |
| agev[m][i]=1; |
agev[m][i]=1; |
|
Line 4733 Title=%s <br>Datafile=%s Firstpass=%d La
|
Line 4929 Title=%s <br>Datafile=%s Firstpass=%d La
|
| |
|
| /* Calculates basic frequencies. Computes observed prevalence at single age |
/* Calculates basic frequencies. Computes observed prevalence at single age |
| and prints on file fileres'p'. */ |
and prints on file fileres'p'. */ |
| freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint); |
freqsummary(fileres, agemin, agemax, s, agev, nlstate, imx,Tvaraff,nbcode, ncodemax,mint,anint,strstart); |
| |
|
| fprintf(fichtm,"\n"); |
fprintf(fichtm,"\n"); |
| fprintf(fichtm,"<br>Total number of observations=%d <br>\n\ |
fprintf(fichtm,"<br>Total number of observations=%d <br>\n\ |
|
Line 4752 Interval (in months) between two waves:
|
Line 4948 Interval (in months) between two waves:
|
| p=param[1][1]; /* *(*(*(param +1)+1)+0) */ |
p=param[1][1]; /* *(*(*(param +1)+1)+0) */ |
| |
|
| globpr=0; /* To get the number ipmx of contributions and the sum of weights*/ |
globpr=0; /* To get the number ipmx of contributions and the sum of weights*/ |
| |
|
| if (mle==-3){ |
if (mle==-3){ |
| ximort=matrix(1,NDIM,1,NDIM); |
ximort=matrix(1,NDIM,1,NDIM); |
| cens=ivector(1,n); |
cens=ivector(1,n); |
|
Line 4761 Interval (in months) between two waves:
|
Line 4958 Interval (in months) between two waves:
|
| |
|
| for (i=1; i<=imx; i++){ |
for (i=1; i<=imx; i++){ |
| dcwave[i]=-1; |
dcwave[i]=-1; |
| for (j=1; j<=lastpass; j++) |
for (m=firstpass; m<=lastpass; m++) |
| if (s[j][i]>nlstate) { |
if (s[m][i]>nlstate) { |
| dcwave[i]=j; |
dcwave[i]=m; |
| /* printf("i=%d j=%d s=%d dcwave=%d\n",i,j, s[j][i],dcwave[i]);*/ |
/* printf("i=%d j=%d s=%d dcwave=%d\n",i,j, s[j][i],dcwave[i]);*/ |
| break; |
break; |
| } |
} |
|
Line 4772 Interval (in months) between two waves:
|
Line 4969 Interval (in months) between two waves:
|
| for (i=1; i<=imx; i++) { |
for (i=1; i<=imx; i++) { |
| if (wav[i]>0){ |
if (wav[i]>0){ |
| ageexmed[i]=agev[mw[1][i]][i]; |
ageexmed[i]=agev[mw[1][i]][i]; |
| j=wav[i];agecens[i]=1.; |
j=wav[i]; |
| if (ageexmed[i]>1 & wav[i]>0) agecens[i]=agev[mw[j][i]][i]; |
agecens[i]=1.; |
| cens[i]=1; |
|
| |
if (ageexmed[i]> 1 && wav[i] > 0){ |
| if (ageexmed[i]<1) cens[i]=-1; |
agecens[i]=agev[mw[j][i]][i]; |
| if (agedc[i]< AGESUP & agedc[i]>1 & dcwave[i]>firstpass & dcwave[i]<=lastpass) cens[i]=0 ; |
cens[i]= 1; |
| |
}else if (ageexmed[i]< 1) |
| |
cens[i]= -1; |
| |
if (agedc[i]< AGESUP && agedc[i]>1 && dcwave[i]>firstpass && dcwave[i]<=lastpass) |
| |
cens[i]=0 ; |
| } |
} |
| else cens[i]=-1; |
else cens[i]=-1; |
| } |
} |
|
Line 4786 Interval (in months) between two waves:
|
Line 4987 Interval (in months) between two waves:
|
| for (j=1;j<=NDIM;j++) |
for (j=1;j<=NDIM;j++) |
| ximort[i][j]=(i == j ? 1.0 : 0.0); |
ximort[i][j]=(i == j ? 1.0 : 0.0); |
| } |
} |
| |
|
| p[1]=0.1; p[2]=0.1; |
p[1]=0.0268; p[NDIM]=0.083; |
| /*printf("%lf %lf", p[1], p[2]);*/ |
/*printf("%lf %lf", p[1], p[2]);*/ |
| |
|
| |
|
| printf("Powell\n"); fprintf(ficlog,"Powell\n"); |
printf("Powell\n"); fprintf(ficlog,"Powell\n"); |
| strcpy(filerespow,"pow-mort"); |
strcpy(filerespow,"pow-mort"); |
| strcat(filerespow,fileres); |
strcat(filerespow,fileres); |
| if((ficrespow=fopen(filerespow,"w"))==NULL) { |
if((ficrespow=fopen(filerespow,"w"))==NULL) { |
| printf("Problem with resultfile: %s\n", filerespow); |
printf("Problem with resultfile: %s\n", filerespow); |
| fprintf(ficlog,"Problem with resultfile: %s\n", filerespow); |
fprintf(ficlog,"Problem with resultfile: %s\n", filerespow); |
| } |
} |
| fprintf(ficrespow,"# Powell\n# iter -2*LL"); |
fprintf(ficrespow,"# Powell\n# iter -2*LL"); |
| /* 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++) |
| if(j!=i)fprintf(ficrespow," p%1d%1d",i,j); |
if(j!=i)fprintf(ficrespow," p%1d%1d",i,j); |
| */ |
*/ |
| fprintf(ficrespow,"\n"); |
fprintf(ficrespow,"\n"); |
| |
|
| powell(p,ximort,NDIM,ftol,&iter,&fret,gompertz); |
powell(p,ximort,NDIM,ftol,&iter,&fret,gompertz); |
| fclose(ficrespow); |
fclose(ficrespow); |
| |
|
| hesscov(matcov, p, NDIM,delti, 1e-4, gompertz); |
hesscov(matcov, p, NDIM, delti, 1e-4, gompertz); |
| |
|
| for(i=1; i <=NDIM; i++) |
for(i=1; i <=NDIM; i++) |
| for(j=i+1;j<=NDIM;j++) |
for(j=i+1;j<=NDIM;j++) |
|
Line 4825 Interval (in months) between two waves:
|
Line 5026 Interval (in months) between two waves:
|
| printf("iter=%d MLE=%f Eq=%lf*exp(%lf*(age-%d))\n",iter,-gompertz(p),p[1],p[2],agegomp); |
printf("iter=%d MLE=%f Eq=%lf*exp(%lf*(age-%d))\n",iter,-gompertz(p),p[1],p[2],agegomp); |
| for (i=1;i<=NDIM;i++) |
for (i=1;i<=NDIM;i++) |
| printf("%f [%f ; %f]\n",p[i],p[i]-2*sqrt(matcov[i][i]),p[i]+2*sqrt(matcov[i][i])); |
printf("%f [%f ; %f]\n",p[i],p[i]-2*sqrt(matcov[i][i]),p[i]+2*sqrt(matcov[i][i])); |
| |
|
| |
lsurv=vector(1,AGESUP); |
| |
lpop=vector(1,AGESUP); |
| |
tpop=vector(1,AGESUP); |
| |
lsurv[agegomp]=100000; |
| |
|
| |
for (k=agegomp;k<=AGESUP;k++) { |
| |
agemortsup=k; |
| |
if (p[1]*exp(p[2]*(k-agegomp))>1) break; |
| |
} |
| |
|
| |
for (k=agegomp;k<agemortsup;k++) |
| |
lsurv[k+1]=lsurv[k]-lsurv[k]*(p[1]*exp(p[2]*(k-agegomp))); |
| |
|
| |
for (k=agegomp;k<agemortsup;k++){ |
| |
lpop[k]=(lsurv[k]+lsurv[k+1])/2.; |
| |
sumlpop=sumlpop+lpop[k]; |
| |
} |
| |
|
| |
tpop[agegomp]=sumlpop; |
| |
for (k=agegomp;k<(agemortsup-3);k++){ |
| |
/* tpop[k+1]=2;*/ |
| |
tpop[k+1]=tpop[k]-lpop[k]; |
| |
} |
| |
|
| |
|
| |
printf("\nAge lx qx dx Lx Tx e(x)\n"); |
| |
for (k=agegomp;k<(agemortsup-2);k++) |
| |
printf("%d %.0lf %lf %.0lf %.0lf %.0lf %lf\n",k,lsurv[k],p[1]*exp(p[2]*(k-agegomp)),(p[1]*exp(p[2]*(k-agegomp)))*lsurv[k],lpop[k],tpop[k],tpop[k]/lsurv[k]); |
| |
|
| |
|
| replace_back_to_slash(pathc,path); /* Even gnuplot wants a / */ |
replace_back_to_slash(pathc,path); /* Even gnuplot wants a / */ |
| printinggnuplotmort(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p); |
printinggnuplotmort(fileres, optionfilefiname,ageminpar,agemaxpar,fage, pathc,p); |
| |
|
| printinghtmlmort(fileres,title,datafile, firstpass, lastpass, \ |
printinghtmlmort(fileres,title,datafile, firstpass, lastpass, \ |
| stepm, weightopt,\ |
stepm, weightopt,\ |
| model,imx,p,matcov); |
model,imx,p,matcov,agemortsup); |
| |
|
| |
free_vector(lsurv,1,AGESUP); |
| |
free_vector(lpop,1,AGESUP); |
| |
free_vector(tpop,1,AGESUP); |
| } /* Endof if mle==-3 */ |
} /* Endof if mle==-3 */ |
| |
|
| else{ /* For mle >=1 */ |
else{ /* For mle >=1 */ |
| |
|
| likelione(ficres, p, npar, nlstate, &globpr, &ipmx, &sw, &fretone, funcone); /* Prints the contributions to the likelihood */ |
likelione(ficres, p, npar, nlstate, &globpr, &ipmx, &sw, &fretone, funcone); /* Prints the contributions to the likelihood */ |
|
Line 5090 Interval (in months) between two waves:
|
Line 5326 Interval (in months) between two waves:
|
| } |
} |
| printf("Computing stable prevalence: result on file '%s' \n", filerespl); |
printf("Computing stable prevalence: result on file '%s' \n", filerespl); |
| fprintf(ficlog,"Computing stable prevalence: result on file '%s' \n", filerespl); |
fprintf(ficlog,"Computing stable prevalence: result on file '%s' \n", filerespl); |
| |
fprintf(ficrespl, "#Local time at start: %s", strstart); |
| fprintf(ficrespl,"#Stable prevalence \n"); |
fprintf(ficrespl,"#Stable prevalence \n"); |
| fprintf(ficrespl,"#Age "); |
fprintf(ficrespl,"#Age "); |
| for(i=1; i<=nlstate;i++) fprintf(ficrespl,"%d-%d ",i,i); |
for(i=1; i<=nlstate;i++) fprintf(ficrespl,"%d-%d ",i,i); |
|
Line 5150 Interval (in months) between two waves:
|
Line 5387 Interval (in months) between two waves:
|
| hstepm=hstepm/stepm; /* Typically 2 years, = 2/6 months = 4 */ |
hstepm=hstepm/stepm; /* Typically 2 years, = 2/6 months = 4 */ |
| |
|
| /* hstepm=1; aff par mois*/ |
/* hstepm=1; aff par mois*/ |
| |
fprintf(ficrespij, "#Local time at start: %s", strstart); |
| fprintf(ficrespij,"#****** h Pij x Probability to be in state j at age x+h being in i at x "); |
fprintf(ficrespij,"#****** h Pij x Probability to be in state j at age x+h being in i at x "); |
| for(cptcov=1,k=0;cptcov<=i1;cptcov++){ |
for(cptcov=1,k=0;cptcov<=i1;cptcov++){ |
| for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){ |
for(cptcod=1;cptcod<=ncodemax[cptcov];cptcod++){ |
|
Line 5187 Interval (in months) between two waves:
|
Line 5424 Interval (in months) between two waves:
|
| } |
} |
| } |
} |
| |
|
| varprob(optionfilefiname, matcov, p, delti, nlstate, bage, fage,k,Tvar,nbcode, ncodemax); |
varprob(optionfilefiname, matcov, p, delti, nlstate, bage, fage,k,Tvar,nbcode, ncodemax,strstart); |
| |
|
| fclose(ficrespij); |
fclose(ficrespij); |
| |
|
|
Line 5276 Interval (in months) between two waves:
|
Line 5513 Interval (in months) between two waves:
|
| |
|
| eij=ma3x(1,nlstate,1,nlstate,(int) bage, (int) fage); |
eij=ma3x(1,nlstate,1,nlstate,(int) bage, (int) fage); |
| oldm=oldms;savm=savms; |
oldm=oldms;savm=savms; |
| evsij(fileres, eij, p, nlstate, stepm, (int) bage, (int)fage, oldm, savm, k, estepm, delti, matcov); |
evsij(fileres, eij, p, nlstate, stepm, (int) bage, (int)fage, oldm, savm, k, estepm, delti, matcov, strstart); |
| |
|
| vareij=ma3x(1,nlstate,1,nlstate,(int) bage, (int) fage); |
vareij=ma3x(1,nlstate,1,nlstate,(int) bage, (int) fage); |
| oldm=oldms;savm=savms; |
oldm=oldms;savm=savms; |
| varevsij(optionfilefiname, vareij, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl,k, estepm, cptcov,cptcod,0, mobilav); |
varevsij(optionfilefiname, vareij, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl,k, estepm, cptcov,cptcod,0, mobilav, strstart); |
| if(popbased==1){ |
if(popbased==1){ |
| varevsij(optionfilefiname, vareij, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl,k, estepm, cptcov,cptcod,popbased,mobilav); |
varevsij(optionfilefiname, vareij, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl,k, estepm, cptcov,cptcod,popbased,mobilav, strstart); |
| } |
} |
| |
|
| |
fprintf(ficrest, "#Local time at start: %s", strstart); |
| fprintf(ficrest,"#Total LEs with variances: e.. (std) "); |
fprintf(ficrest,"#Total LEs with variances: e.. (std) "); |
| for (i=1;i<=nlstate;i++) fprintf(ficrest,"e.%d (std) ",i); |
for (i=1;i<=nlstate;i++) fprintf(ficrest,"e.%d (std) ",i); |
| fprintf(ficrest,"\n"); |
fprintf(ficrest,"\n"); |
|
Line 5358 Interval (in months) between two waves:
|
Line 5595 Interval (in months) between two waves:
|
| |
|
| varpl=matrix(1,nlstate,(int) bage, (int) fage); |
varpl=matrix(1,nlstate,(int) bage, (int) fage); |
| oldm=oldms;savm=savms; |
oldm=oldms;savm=savms; |
| varprevlim(fileres, varpl, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl,k); |
varprevlim(fileres, varpl, matcov, p, delti, nlstate, stepm, (int) bage, (int) fage, oldm, savm, prlim, ftolpl,k,strstart); |
| free_matrix(varpl,1,nlstate,(int) bage, (int)fage); |
free_matrix(varpl,1,nlstate,(int) bage, (int)fage); |
| } |
} |
| } |
} |
|
Line 5407 Interval (in months) between two waves:
|
Line 5644 Interval (in months) between two waves:
|
| tm = *localtime(&end_time.tv_sec); |
tm = *localtime(&end_time.tv_sec); |
| tmg = *gmtime(&end_time.tv_sec); |
tmg = *gmtime(&end_time.tv_sec); |
| strcpy(strtend,asctime(&tm)); |
strcpy(strtend,asctime(&tm)); |
| printf("Local time at start %s\nLocaltime at end %s",strstart, strtend); |
printf("Local time at start %s\nLocal time at end %s",strstart, strtend); |
| fprintf(ficlog,"Local time at start %s\nLocal time at end %s\n",strstart, strtend); |
fprintf(ficlog,"Local time at start %s\nLocal time at end %s\n",strstart, strtend); |
| printf("Total time used %s\n", asc_diff_time(end_time.tv_sec -start_time.tv_sec,tmpout)); |
printf("Total time used %s\n", asc_diff_time(end_time.tv_sec -start_time.tv_sec,tmpout)); |
| |
|
|
Line 5424 Interval (in months) between two waves:
|
Line 5661 Interval (in months) between two waves:
|
| /*------ End -----------*/ |
/*------ End -----------*/ |
| |
|
| chdir(path); |
chdir(path); |
| strcpy(plotcmd,"\""); |
/*strcat(plotcmd,CHARSEPARATOR);*/ |
| strcat(plotcmd,pathimach); |
sprintf(plotcmd,"gnuplot"); |
| strcat(plotcmd,GNUPLOTPROGRAM); |
#ifndef UNIX |
| strcat(plotcmd,"\""); |
sprintf(plotcmd,"\"%swgnuplot.exe\"",pathimach); |
| strcat(plotcmd," "); |
#endif |
| strcat(plotcmd,optionfilegnuplot); |
if(!stat(plotcmd,&info)){ |
| printf("Starting graphs with: %s",plotcmd);fflush(stdout); |
printf("Error gnuplot program not found: %s\n",plotcmd);fflush(stdout); |
| |
if(!stat(getenv("GNUPLOTBIN"),&info)){ |
| |
printf("Error gnuplot program not found: %s Environment GNUPLOTBIN not set.\n",plotcmd);fflush(stdout); |
| |
}else |
| |
strcpy(pplotcmd,plotcmd); |
| |
#ifdef UNIX |
| |
strcpy(plotcmd,GNUPLOTPROGRAM); |
| |
if(!stat(plotcmd,&info)){ |
| |
printf("Error gnuplot program not found: %s\n",plotcmd);fflush(stdout); |
| |
}else |
| |
strcpy(pplotcmd,plotcmd); |
| |
#endif |
| |
}else |
| |
strcpy(pplotcmd,plotcmd); |
| |
|
| |
sprintf(plotcmd,"%s %s",pplotcmd, optionfilegnuplot); |
| |
printf("Starting graphs with: %s\n",plotcmd);fflush(stdout); |
| |
|
| if((outcmd=system(plotcmd)) != 0){ |
if((outcmd=system(plotcmd)) != 0){ |
| printf(" Problem with gnuplot\n"); |
printf("\n Problem with gnuplot\n"); |
| } |
} |
| printf(" Wait..."); |
printf(" Wait..."); |
| while (z[0] != 'q') { |
while (z[0] != 'q') { |
![[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