Dear statalists,

I tried to conduct a Monte Carlo simulation to show that under conditional heteroscedasticity, gmm estimation of the probit model is consistent, but my simulation result doesn't support this (shown in the graph). The heteroscedasticity is in a multiplicative form described in Harvey (1976), and it doesn't seem I got the moment condition wrong. I also simulate the estimation from -probit- and -hetprobit-, which, as expected, are inconsistent and consistent, respectively. What might be the problem here? Thanks!

Array

Code:
clear all
mkdir d:/temp
cd d:/temp
cap log close
log using probit_hetroscedasticity.log, replace

* generate simulation data and estimate the model using ols
cap program drop ph
program define ph, rclass
    drop _all
    syntax [, obs(integer 1)]
    set obs `obs'
    
    * DGP with conditional heteroscedasticity
    gen double x = rnormal()
    gen double e = rnormal() 
    gen double v2 = exp(0.4*x) * e
    gen double y_star = 1 + 2*x + v2  
    gen y = y_star > 0
    
    probit y x, iterate(50)
    scalar b_homo = _b[x]
    
    hetprobit y x, het(x)
    scalar b_het = _b[x]
    
    gmm (y - normal({b0} + {b1}*x)), instruments(x)
    mat beta = e(b)
    scalar b_gmm = beta[1,2]
end

simulate b_homo = b_homo b_het = b_het b_gmm = b_gmm, ///
         reps(500) saving(d1, replace) : ph, obs(100)

simulate b_homo2 = b_homo b_het2 = b_het b_gmm2 = b_gmm, ///
         reps(500) saving(d2, replace) : ph, obs(1000)

* merge all estimated parameters
use d1, clear
merge 1:1 _n using d2
drop _merge
rm d1.dta 
rm d2.dta

* Generate density graph of the estimated parameters 
kdensity b_homo, gen(b_homo_x b_homo_d) nograph n(500)
kdensity b_het, gen(b_het_x b_het_d) nograph n(500)
kdensity b_gmm, gen(b_gmm_x b_gmm_d) nograph n(500)

kdensity b_homo2, gen(b_homo2_x b_homo2_d) nograph n(500)
kdensity b_het2, gen(b_het2_x b_het2_d) nograph n(500)
kdensity b_gmm2, gen(b_gmm2_x b_gmm2_d) nograph n(500)


label var b_homo "N=100, homoscedasticity"
label var b_het "N=100, heteroscedasticity"
label var b_gmm "N=100, GMM"

label var b_homo2 "N=1000, homoscedasticity"
label var b_het2 "N=1000, heteroscedasticity"
label var b_gmm2 "N=1000, GMM"

save data, replace

* generate pdf
graph twoway (line b_homo_d b_homo_x) (line b_het_d b_het_x) (line b_gmm_d b_gmm_x) ///
             (line b_homo2_d b_homo2_x) (line b_het2_d b_het2_x) (line b_gmm2_d b_gmm2_x), xline(2, lcolor(red))
graph export probit_estimation_heteroscedasticity.png, replace

log close