question regarding ML estimation

[Sungmin Cheu]

Dear Statalists,

I am trying to bivariate regression and now working on ml program.

Shortened version of the equations that I try to estmate are as below.

Equatian1: log_business_income = beta0 + beta1 * subsidy + beta2 * age
Equation2: log_labor_income = beta0 + beta1 * subsidy + beta2 * age


I set up likelihood as shown in the below code.
I have three situations where 1) an individual only participates in business, 2) only in labor, and 3) both in business and labor.
However, ml check keeps giving me 'could not find feasible values' message and I am without any clue about why so.
I assume that I correctly set up the likelihood function.

Please enlighten me if there is anything wrong with the code below...

* theta3: std. deviation of log_business_inc / theta4: std. deviation of log_labor_inc / theta5: correlation between log_business_inc and log_labor_inc
************************************************** ************
program part2
version 14.0
args todo b lnf
tempvar theta1 theta2 theta3 theta4 theta5 lnfj
mleval `theta1'=`b', eq(1)
mleval `theta2'=`b', eq(2)
mleval `theta3'=`b', eq(3)
mleval `theta4'=`b', eq(4)
mleval `theta5'=`b', eq(5)
gen double `lnfj'=ln(1/sqrt(2*`theta3'^2)) - (1/2)*((($ML_y1 - `theta1') / `theta3')^2) if $ML_y1>0 & $ML_y2==.
replace `lnfj'=ln(1/sqrt(2*`theta4'^2)) - (1/2)*((($ML_y2 - `theta2') / `theta4')^2) if $ML_y1==. & $ML_y2>0
replace `lnfj'=ln(1/sqrt(`theta3'^2*`theta4'^2*(1 - `theta5'^2))) - .5*(1/(1-`theta5'^2))*((($ML_y1 - `theta1')/`theta3')^2 + (($ML_y2 - `theta2')/`theta4')^2 - 2*`theta5'*(($ML_y1 - `theta1')/`theta3')*(($ML_y2 - `theta2')/`theta4')) if $ML_y1>0 & $ML_y2==.
mlsum `lnf'=`lnfj'
end

ml model d0 part2 (xb: log_business_inc = g_sub age) (xc: log_labor_inc = g_sub age ) () () ()
ml check
ml search
ml maximize
************************************************** ************

Thanking you in advance.
Cheers.

[Matthew J. Baker]

Sungmin Cheu --

The only things that I see that might be problems are 1) you use `theta3'^2 as a standard deviation parameter, which is fine, but I have found that sometimes it is better to use an exponentiated version to avoid the possibility of negative numbers - i.e., use exp('theta3') instead.

Also, I'm not sure if this is intentional, but it looks like you use the condition:

 if $ML_y1>0 & $ML_y2==.

twice when generating lnfj - both in the initial generate statement, and then in the last replace statement.

Just my two cents --

Matthew J. Baker