A Grid Search for maximum likelihood estimations

[Jinrui Pan]

Dear All,

I am fairly new to Stata Programming. I am wondering if anyone can help me code a grid search of ml estimation. I've written the program and did the estimation, however, I've encountered some problems when adding treatments - I'd like to obtain estimations of the parameters in different treatment groups. Therefore, I want to try a grid search for appropriate parameter values but not sure where to start... I did a fair amount of search online, and haven't found any solutions, especially for 3 parameters estimation.

Here is my program code:
program define reinforcement
args logl LGphi LGeps LNlam

quietly{
// Transform parameters to respect bounds
tempvar lam phi eps
gen double `lam' = exp(`LNlam')
gen double `phi' = logistic(`LGphi')
gen double `eps' = logistic(`LGeps')

// Initialise values of attraction to each strategy
replace A1 = .
replace A2 = .
replace A3 = .
replace A4 = .
replace A5 = .
replace A6 = .
replace A12 = .

replace cf1 = (1-`eps')*r1 if s == 1
replace cf1 = (`eps'/6)*r1 if s != 1
replace cf2 = (1-`eps')*r2 if s == 2
replace cf2 = (`eps'/6)*r2 if s != 2
replace cf3 = (1-`eps')*r3 if s == 3
replace cf3 = (`eps'/6)*r3 if s != 3
replace cf4 = (1-`eps')*r4 if s == 4
replace cf4 = (`eps'/6)*r4 if s != 4
replace cf5 = (1-`eps')*r5 if s == 5
replace cf5 = (`eps'/6)*r5 if s != 5
replace cf6 = (1-`eps')*r6 if s == 6
replace cf6 = (`eps'/6)*r6 if s != 6
replace cf12 = (1-`eps')*r12 if s == 12
replace cf12 = (`eps'/6)*r12 if s != 12

// Fill in the value of attraction for first obs on each subject
by id: replace A1 = (1-`phi')*0.1613 + cf1 if _n == 1
by id: replace A1 = (1-`phi')*0.1613 + cf1 if _n == 1
by id: replace A2 = (1-`phi')*0.1548 + cf2 if _n == 1
by id: replace A2 = (1-`phi')*0.1548 + cf2 if _n == 1
by id: replace A3 = (1-`phi')*0.1484 + cf3 if _n == 1
by id: replace A3 = (1-`phi')*0.1484 + cf3 if _n == 1
by id: replace A4 = (1-`phi')*0.1677 + cf4 if _n == 1
by id: replace A4 = (1-`phi')*0.1677 + cf4 if _n == 1
by id: replace A5 = (1-`phi')*0.0387 + cf5 if _n == 1
by id: replace A5 = (1-`phi')*0.0387 + cf5 if _n == 1
by id: replace A6 = (1-`phi')*0.1871 + cf6 if _n == 1
by id: replace A6 = (1-`phi')*0.1871 + cf6 if _n == 1
by id: replace A12 = (1-`phi')*0.1419 + cf12 if _n == 1
by id: replace A12 = (1-`phi')*0.1419 + cf12 if _n == 1

// Update the value of attraction for subsequent observations
by id: replace A1 = (1-`phi')*A1[_n-1] + cf1 if A1 == . & _n > 1
by id: replace A1 = (1-`phi')*A1[_n-1] + cf1 if A1 == . & _n > 1
by id: replace A2 = (1-`phi')*A2[_n-1] + cf2 if A2 == . & _n > 1
by id: replace A2 = (1-`phi')*A2[_n-1] + cf2 if A2 == . & _n > 1
by id: replace A3 = (1-`phi')*A3[_n-1] + cf3 if A3 == . & _n > 1
by id: replace A3 = (1-`phi')*A3[_n-1] + cf3 if A3 == . & _n > 1
by id: replace A4 = (1-`phi')*A4[_n-1] + cf4 if A4 == . & _n > 1
by id: replace A4 = (1-`phi')*A4[_n-1] + cf4 if A4 == . & _n > 1
by id: replace A5 = (1-`phi')*A5[_n-1] + cf5 if A5 == . & _n > 1
by id: replace A5 = (1-`phi')*A5[_n-1] + cf5 if A5 == . & _n > 1
by id: replace A6 = (1-`phi')*A6[_n-1] + cf6 if A6 == . & _n > 1
by id: replace A6 = (1-`phi')*A6[_n-1] + cf6 if A6 == . & _n > 1
by id: replace A12 = (1-`phi')*A12[_n-1] + cf12 if A12 == . & _n > 1
by id: replace A12 = (1-`phi')*A12[_n-1] + cf12 if A12 == . & _n > 1

*noi sum A2-A5
// Initialise the probability of choosing strategy 4
replace p1 = .
replace p2 = .
replace p3 = .
replace p4 = .
replace p5 = .
replace p6 = .

// Predict probability of choosing strategy 4 under logit and given parameter estimates
by id: replace p1 = exp(`lam'*.1050622)/ ///
(exp(`lam'*.1050622)+exp(`lam'*.0639338)+exp(`lam' *.0217152)+ ///
exp(`lam'*.1439694)+exp(`lam'*-1.322084)+exp(`lam'*.2534271)+ ///
exp(`lam'* -.0230691)) if _n == 1
by id: replace p2 = exp(`lam'*.0639338)/ ///
(exp(`lam'*.1050622)+exp(`lam'*.0639338)+exp(`lam' *.0217152)+ ///
exp(`lam'*.1439694)+exp(`lam'*-1.322084)+exp(`lam'*.2534271)+ ///
exp(`lam'* -.0230691)) if _n == 1
by id: replace p3 = exp(`lam'*.0217152)/ ///
(exp(`lam'*.1050622)+exp(`lam'*.0639338)+exp(`lam' *.0217152)+ ///
exp(`lam'*.1439694)+exp(`lam'*-1.322084)+exp(`lam'*.2534271)+ ///
exp(`lam'* -.0230691)) if _n == 1
by id: replace p4 = exp(`lam'*.1439694)/ ///
(exp(`lam'*.1050622)+exp(`lam'*.0639338)+exp(`lam' *.0217152)+ ///
exp(`lam'*.1439694)+exp(`lam'*-1.322084)+exp(`lam'*.2534271)+ ///
exp(`lam'* -.0230691)) if _n == 1
by id: replace p5 = exp(`lam'*-1.322084)/ ///
(exp(`lam'*.1050622)+exp(`lam'*.0639338)+exp(`lam' *.0217152)+ ///
exp(`lam'*.1439694)+exp(`lam'*-1.322084)+exp(`lam'*.2534271)+ ///
exp(`lam'* -.0230691)) if _n == 1
by id: replace p6 = exp(`lam'*.2534271)/ ///
(exp(`lam'*.1050622)+exp(`lam'*.0639338)+exp(`lam' *.0217152)+ ///
exp(`lam'*.1439694)+exp(`lam'*-1.322084)+exp(`lam'*.2534271)+ ///
exp(`lam'* -.0230691)) if _n == 1

by id: replace p1 = exp(`lam'*A1)/ ///
(exp(`lam'*A1)+exp(`lam'*A2)+exp(`lam'*A3)+exp(`la m'*A4)+ ///
exp(`lam'*A5)+exp(`lam'*A6)+exp(`lam'*A12)) if p1 == . & _n > 1
by id: replace p2 = exp(`lam'*A2)/ ///
(exp(`lam'*A1)+exp(`lam'*A2)+exp(`lam'*A3)+exp(`la m'*A4)+ ///
exp(`lam'*A5)+exp(`lam'*A6)+exp(`lam'*A12)) if p2 == . & _n > 1
by id: replace p3 = exp(`lam'*A3)/ ///
(exp(`lam'*A1)+exp(`lam'*A2)+exp(`lam'*A3)+exp(`la m'*A4)+ ///
exp(`lam'*A5)+exp(`lam'*A6)+exp(`lam'*A12)) if p3 == . & _n > 1
by id: replace p4 = exp(`lam'*A4)/ ///
(exp(`lam'*A1)+exp(`lam'*A2)+exp(`lam'*A3)+exp(`la m'*A4)+ ///
exp(`lam'*A5)+exp(`lam'*A6)+exp(`lam'*A12)) if p4 == . & _n > 1
by id: replace p5 = exp(`lam'*A5)/ ///
(exp(`lam'*A1)+exp(`lam'*A2)+exp(`lam'*A3)+exp(`la m'*A4)+ ///
exp(`lam'*A5)+exp(`lam'*A6)+exp(`lam'*A12)) if p5 == . & _n > 1
by id: replace p6 = exp(`lam'*A6)/ ///
(exp(`lam'*A1)+exp(`lam'*A2)+exp(`lam'*A3)+exp(`la m'*A4)+ ///
exp(`lam'*A5)+exp(`lam'*A6)+exp(`lam'*A12)) if p6 == . & _n > 1

replace `logl'=ln(p1)*s1+ln(p2)*s2+ln(p3)*s3+ln(p4)*s4+ln( p5)*s5+ln(p6)*s6+ln(`p12')*(1-s1-s2-s3-s4-s5-s6)
*/
replace `logl'=ln(p1)*s1+ln(p2)*s2+ln(p3)*s3+ln(p4)*s4+ln( p5)*s5+ln(p6)*s6+ln(1-p1-p2-p3-p4-p5-p6)*(1-s1-s2-s3-s4-s5-s6)
}
end

Any help would be massively appreciated!!!