From 1fd29804047c48aa064111c7a30d78f656856ce7 Mon Sep 17 00:00:00 2001 From: NBrouard Date: Thu, 26 Dec 2024 16:33:37 +0100 Subject: [PATCH] * src/imach.c: improving warning on Hessian inversion and nan --- src/imach.c | 33 +++++++++++++++++---------------- 1 file changed, 17 insertions(+), 16 deletions(-) diff --git a/src/imach.c b/src/imach.c index ef8910f..d70ee5f 100644 --- a/src/imach.c +++ b/src/imach.c @@ -1427,7 +1427,7 @@ double gnuplotversion=GNUPLOTVERSION; /* $State$ */ #include "version.h" char version[]=__IMACH_VERSION__; -char copyright[]="November 2024,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2020, Nihon University 2021-202, INED 2000-2024"; +char copyright[]="December 2024,INED-EUROREVES-Institut de longevite-Japan Society for the Promotion of Science (Grant-in-Aid for Scientific Research 25293121), Intel Software 2015-2020, Nihon University 2021-202, INED 2000-2024"; char fullversion[]="$Revision$ $Date$"; char strstart[80]; char optionfilext[10], optionfilefiname[FILENAMELENGTH]; @@ -6785,6 +6785,18 @@ void hesscov(double **matcov, double **hess, double p[], int npar, double delti[ printf("\n"); fprintf(ficlog,"\n"); + + printf("\n#Hessian matrix#\n"); + fprintf(ficlog,"\n#Hessian matrix#\n"); + for (i=1;i<=npar;i++) { + for (j=1;j<=npar;j++) { + printf("%.6e ",hess[i][j]); + fprintf(ficlog,"%.6e ",hess[i][j]); + } + printf("\n"); + fprintf(ficlog,"\n"); + } + printf("\nInverting the hessian to get the covariance matrix. Wait...\n"); fprintf(ficlog,"\nInverting the hessian to get the covariance matrix. Wait...\n"); @@ -6805,17 +6817,6 @@ void hesscov(double **matcov, double **hess, double p[], int npar, double delti[ } } - printf("\n#Hessian matrix#\n"); - fprintf(ficlog,"\n#Hessian matrix#\n"); - for (i=1;i<=npar;i++) { - for (j=1;j<=npar;j++) { - printf("%.6e ",hess[i][j]); - fprintf(ficlog,"%.6e ",hess[i][j]); - } - printf("\n"); - fprintf(ficlog,"\n"); - } - /* printf("\n#Covariance matrix#\n"); */ /* fprintf(ficlog,"\n#Covariance matrix#\n"); */ /* for (i=1;i<=npar;i++) { */ @@ -6829,7 +6830,7 @@ void hesscov(double **matcov, double **hess, double p[], int npar, double delti[ /* Recompute Inverse */ /* for (i=1;i<=npar;i++) */ - /* for (j=1;j<=npar;j++) a[i][j]=matcov[i][j]; */ + /* for (j=1;j<=npar;j+<+) a[i][j]=matcov[i][j]; */ /* ludcmp(a,npar,indx,&pd); */ /* printf("\n#Hessian matrix recomputed#\n"); */ @@ -16244,9 +16245,9 @@ Please run with mle=-1 to get a correct covariance matrix.\n",ageminpar,agemaxpa if(mle != 0){ /* Computing hessian and covariance matrix only at a peak of the Likelihood, that is after optimization */ - printf("Parameters and 95%% confidence intervals\n W is simply the result of the division of the parameter by the square root of covariance of the parameter.\n And Wald-based confidence intervals plus and minus 1.96 * W .\n But be careful that parameters are highly correlated because incidence of disability is highly correlated to incidence of recovery.\n It might be better to visualize the covariance matrix. See the page 'Matrix of variance-covariance of one-step probabilities' and its graphs.\nPlease notice that if the LogLikelihood is not maximised in all directions, the Hessian should not be inverted; but sometimes it can be successfully inverted but producing wrong results like negative variances. In those case the square roots of a negative number produce 'nan' which means 'not a number'.\n"); - fprintf(ficlog, "Parameters, Wald tests and Wald-based confidence intervals\n W is simply the result of the division of the parameter by the square root of covariance of the parameter.\n And Wald-based confidence intervals plus and minus 1.96 * W \n It might be better to visualize the covariance matrix. See the page 'Matrix of variance-covariance of one-step probabilities' and its graphs.\nPlease notice that if the LogLikelihood is not maximised in all directions, the Hessian should not be inverted; but sometimes it can be successfully inverted but producing wrong results like negative variances. In those case the square roots of a negative number produce 'nan' which means 'not a number'.\n"); - fprintf(fichtm, "\n

The Wald test results are output only if the maximimzation of the Likelihood is performed (mle=1)\n
Parameters, Wald tests and Wald-based confidence intervals\n
W is simply the result of the division of the parameter by the square root of covariance of the parameter.\n
And Wald-based confidence intervals plus and minus 1.96 * W \n
It might be better to visualize the covariance matrix. See the page 'Matrix of variance-covariance of one-step probabilities and its graphs'.\n
Please notice that if the LogLikelihood is not maximised in all directions, the Hessian should not be inverted; but sometimes it can be successfully inverted but producing wrong results like negative variances. In those case the square roots of a negative number produce 'nan' which means 'not a number'.\n
",optionfilehtmcov); + printf("Parameters and 95%% confidence intervals\n W is simply the result of the division of the parameter by the square root of covariance of the parameter.\n And Wald-based confidence intervals plus and minus 1.96 * W .\n But be careful that parameters are highly correlated because incidence of disability is highly correlated to incidence of recovery.\n It might be better to visualize the covariance matrix. See the page 'Matrix of variance-covariance of one-step probabilities' and its graphs.\nPlease notice that if the LogLikelihood is not maximised in all directions, the Hessian might not be invertable; but sometimes it can be successfully inverted but producing wrong results like negative variances. In those case the square roots of a negative number produce 'nan' which means 'not a number' on Unix or -1.#IND on Windows.\n"); + fprintf(ficlog, "Parameters, Wald tests and Wald-based confidence intervals\n W is simply the result of the division of the parameter by the square root of covariance of the parameter.\n And Wald-based confidence intervals plus and minus 1.96 * W \n It might be better to visualize the covariance matrix. See the page 'Matrix of variance-covariance of one-step probabilities' and its graphs.\nPlease notice that if the LogLikelihood is not maximised in all directions, the Hessian might not be invertable; but sometimes it can be successfully inverted but producing wrong results like negative variances. In those case the square roots of a negative number produce 'nan' which means 'not a number' on Unix or -1.#IND on Windows.\n"); + fprintf(fichtm, "\n

The Wald test results are output only if the maximimzation of the Likelihood is performed (mle=1)\n
Parameters, Wald tests and Wald-based confidence intervals\n
W is simply the result of the division of the parameter by the square root of covariance of the parameter.\n
And Wald-based confidence intervals plus and minus 1.96 * W \n
It might be better to visualize the covariance matrix. See the page 'Matrix of variance-covariance of one-step probabilities and its graphs'.\n
Please notice that if the LogLikelihood is not maximised in all directions, the Hessian might not be invertable; but sometimes it can be successfully inverted but producing wrong results like negative variances. In those case the square roots of a negative number produce 'nan' which means 'not a number' on Unix or -1.#IND on Windows.\n
",optionfilehtmcov); fprintf(fichtm,"\n"); fprintf(fichtm, "\n"); /* ncovmodel=2+nbocc(model,'+')+1; */ -- 2.43.0
Model=1+ age