Ordered Probit with Random Coefficient Using Maximum Simulated Likelihood

[Behram Wali]

Dear All,

I am working on coding a relatively simple ordered probit model with a random parameter/random coefficient for one of the explanatory factor. I am using Maximum Simulated Likelihood to average the integral. At a basic level, I am using draws from normal density using the inverse cumulative for the standard deviation of the random parameter. After running ml check on my code, no errors are returned. However, I am struggling with initial values with the following error returned each time "could not find feasible values"

Just to put this into a context, I have successfully coded a random parameter binary probit/logit with "n" numbers of random parameters using the benchmark example in Cameron and Trivedi (2010). However, am struggling with the ordinal case.

Below is my code:

    
    clear all
    sysuse auto
    recode rep78 (2 = 1) (2 = 1) (5 = 1), gen(nx1)
    set seed 19011992
        *Use this to in turn get draws from the normal distribution
        *Create 50 draws (S=50) from the uniform for 
                *each observation in the data at hand(N=74)

capture drop draws*
forvalues i = 1/50 {
    gen draws`i' = runiform()
}

    capture program drop lfoprobitmsl_bw1
        program lfoprobitmsl_bw1
        version 14.2
        args lnf b1 b2 alpha1 ln_sd                             //b1 is constant term and b2, b3 are slopes, alpha1 is threshold
                                                                                //Dont use 'sd' as it will get lost if sd<0
        tempvar sim_fx sim_avgfx                                 //Get temp holders for the simulated densities
        local y "$ML_y1"
        local sd = exp(`ln_sd')                                    //Tranform back to SD
        qui gen `sim_avgfx' = 0
            set seed 19011992
                forvalues d = 1/50 {
                quietly replace `lnf' = ln(normal(-(`b1' + `b2' + `sd'*invnormal(draws`d')*mpg))) if $ML_y1 == 0
                /*contribution to log-likelihood for outcome y=0 respondents*/
                quietly replace `lnf' = ln(normal(-(`b1' + `b2' + `sd'*invnormal(draws`d')*mpg)+`alpha1') - normal(-(`b1' + `b2' + `sd'*invnormal(draws`d')*mpg))) if $ML_y1 == 1
                /*contribution to log-likelihood for outcome y=1 respondents*/
                quietly replace `lnf' = ln(1-normal(-(`b1' + `b2' + `sd'*invnormal(draws`d')*mpg)+`alpha1')) if $ML_y1 == 2
                /*contribution to log-likelihood for outcome y=2 respondents*/
                qui gen `sim_fx' = `lnf'
                qui replace `sim_avgfx' = `sim_avgfx' + `sim_fx'/50
                drop `sim_fx'
            }
        end
    //Let's now calculate the MSL estimator
    ml model lf lfoprobitmsl_bw1 (b1: nx = ) (b2: mpg, nocons) (alpha1:) (ln_sd:)
    ml init 1 1 1 0, copy
    //ml check
    ml maximize, difficult

To generate the error I am getting, I have used Stata's auto.dta. Just for the sake of this example, I have re-coded one of the variable into three categories so that I can use it as my dependent variable (variable nx) in above code. My apologies for that!

My question is: Am I coding the likelihood part correctly or is the issue really with my poor initial values? Have tried a bunch of initial values though.

Any help will be truly appreciated. Thanks for your time!

Sincerely, Behram