Dear all, I attached a sample script and my data. Using the "omit" command after estimating a ordered probit results in a crash on both Ubuntu (64-bit, current cvs) and Win 7(64-bit, version from last week or so). Hope these information are helpful.

On Mon, 7 Dec 2015, Riccardo (Jack) Lucchetti wrote: On Mon, 7 Dec 2015, Artur T. wrote: > > Dear all, >> I attached a sample script and my data. Using the "omit" command after >> estimating a ordered probit results in a crash on both Ubuntu (64-bit, >> current cvs) and Win 7(64-bit, version from last week or so). >> Hope these information are helpful. >> > > It's a numerical problem. For example, if you multiply your eubor_3m > variable (and its lags) by 10, the problem disappears. The same happens > with ifo_r1_dln. > > It looks as if with the original values some weird condition is triggered > such that the Hessian turns out to be numerically singular. I'll look into > this. > Ok, we've got a rather extraordinary combination of factors here. For a start, your initial model has very few observations compared to the size of the model, PLUS many of the regressors are highly collinear, PLUS the orders of magnitude of the original data don't play very well together. As a consequence, identification in your model is very weak. For example, try the following modification of your script: <hansl> open crash_data.gdt --quiet list LHS = arng_ind_r list ARNG = arng_ind_r_7 arng_ind_r_8 arng_ind_r_9 arng_ind_r_10 list IPIDX = ipidx_d_3 ipidx_d_4 ipidx_d_5 ipidx_d_6 list JOBVAC = jobvac_dln_3 jobvac_dln_4 jobvac_dln_5 jobvac_dln_6 list IFO = ifo_r1_dln_1 ifo_r1_dln_2 ifo_r1_dln_3 ifo_r1_dln_4 list EUBOR = eubor_3m_1 eubor_3m_2 eubor_3m_3 eubor_3m_4 list CDAX = cdax_dln_1 cdax_dln_2 cdax_dln_3 cdax_dln_4 list CRP = crp_sprd_1 crp_sprd_2 crp_sprd_3 crp_sprd_4 list AFTGinl = aftg_inl_dln_3 aftg_inl_dln_4 aftg_inl_dln_5 aftg_inl_dln_6 list AFTGaus = aftg_aus_dln_3 aftg_aus_dln_4 aftg_aus_dln_5 aftg_aus_dln_6 list YC10 = yc10y_1 yc10y_2 yc10y_3 yc10y_4 list A = IPIDX IFO loop foreach i A --quiet $i *= 100 endloop list RHS = const ARNG IPIDX JOBVAC IFO EUBOR CDAX CRP AFTGinl AFTGaus YC10 probit LHS RHS matrix H = invpd($vcv) matrix lambda = eigensym(H) print lambda omit RHS --auto=0.05 </hansl> As you can see, simply by rescaling some of the variables, you get your script to work. However, the Hessian of your base model is _very_ badly conditioned (look at the eigenvalues). This results in a number of subtle numerical problems; the crash happens because, at some point, Cholesky decomposition fails in the gretl_invert_symmetric_matrix() function (it must be said, our criteria are relatively picky compared to other packages). I tried a few strategies for auto-scaling, but none of them were convincing. Anyway, it shouldn't just crash, but stop with an error message. It remains to be seen what the best course of action is under circumstances such as these. ------------------------------------------------------- Riccardo (Jack) Lucchetti Dipartimento di Scienze Economiche e Sociali (DiSES) Università Politecnica delle Marche (formerly known as Università di Ancona) r.lucchetti(a)univpm.it http://www2.econ.univpm.it/servizi/hpp/lucchetti ------------------------------------------------------- _______________________________________________ Gretl-devel mailing list Gretl-devel(a)lists.wfu.edu http://lists.wfu.edu/mailman/listinfo/gretl-devel

On Mon, 7 Dec 2015, Artur T. wrote: > Thanks for your quick responses! > The sample script and datafile is actually just an extract from a > rolling-window exercise, and gretl just crashed at this specific window > involving apparently some numerical issues, as you said, Jack. > > Thanks for the hint with the re-scaling; that will help for the moment but > of course gretl should not crash at any point. Upon further investigation: I'm now inclined to think that the root of the problem lies in the function ordered_hessian() in discrete.c, where we simply compute the Hessian by numerically differentiating the analytical score. I distinctly remember that at the time I thought to myself "this'll do for now", but then I forgot about it; of course, in "hard" cases such as yours, this is a major source of numerical instability. When we find some time to write down a C implementation of the analytical Hessian for the ordered probit model, all the problems should be over.

On Mon, 7 Dec 2015, Allin Cottrell wrote: On Mon, 7 Dec 2015, Riccardo (Jack) Lucchetti wrote: > > On Mon, 7 Dec 2015, Artur T. wrote: >> >> Thanks for your quick responses! >>> The sample script and datafile is actually just an extract from a >>> rolling-window exercise, and gretl just crashed at this specific window >>> involving apparently some numerical issues, as you said, Jack. >>> >>> Thanks for the hint with the re-scaling; that will help for the moment >>> but >>> of course gretl should not crash at any point. >>> >> >> Upon further investigation: I'm now inclined to think that the root of >> the problem lies in the function ordered_hessian() in discrete.c, where we >> simply compute the Hessian by numerically differentiating the analytical >> score. I distinctly remember that at the time I thought to myself "this'll >> do for now", but then I forgot about it; of course, in "hard" cases such as >> yours, this is a major source of numerical instability. When we find some >> time to write down a C implementation of the analytical Hessian for the >> ordered probit model, all the problems should be over. >> > > Ah, good idea. In the meantime I've made sure that when the Hessian can't > be inverted the error propagates correctly, stopping the rot before a crash > occurs. > In the mean^2-time, I've used a marginally more sophisticated way to computed the numerical Hessian for ordered probit and committed it to git. Artur's script shouldn't crash any longer. Please tell me if it does. ------------------------------------------------------- Riccardo (Jack) Lucchetti Dipartimento di Scienze Economiche e Sociali (DiSES) Università Politecnica delle Marche (formerly known as Università di Ancona) r.lucchetti(a)univpm.it http://www2.econ.univpm.it/servizi/hpp/lucchetti ------------------------------------------------------- _______________________________________________ Gretl-devel mailing list Gretl-devel(a)lists.wfu.edu http://lists.wfu.edu/mailman/listinfo/gretl-devel

>> Ah, good idea. In the meantime I've made sure that when the Hessian >> can't be inverted the error propagates correctly, stopping the rot >> before a crash occurs. > > In the mean^2-time, I've used a marginally more sophisticated way to > computed the numerical Hessian for ordered probit and committed it to > git. Artur's script shouldn't crash any longer. Please tell me if it > does. It wouldn't crash after my amendment, but after yours it actually works. Nice job! Allin