Bookmark and Share

Notice: On March 31, it was announced that Statalist is moving from an email list to a forum. The old list will shut down on April 23, and its replacement, statalist.org is already up and running.


[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: st: Right Censored ordered probit with Stata 11


From   Urmi Bhattacharya <ub3@indiana.edu>
To   statalist@hsphsun2.harvard.edu
Subject   Re: st: Right Censored ordered probit with Stata 11
Date   Tue, 27 Sep 2011 00:47:02 +0530

Hi,

I am posting the entire code , the simulated data and the output I get
( my apologies for it being too long).

set trace on

. set traced 1

. *Start cprbt.ado ---
. capture program drop coprbt

. program define coprbt , eclass byable(recall)
  1.         version 11.0
  2. *Replay support.
.         if replay() {
  3.                 if "`e(cmd)'" != "coprbt" {
  4.                         error 301
  5.                 }
  6.                 syntax [, level(passthru)]
  7.         }
  8. *Main program
.         else {
  9.                 syntax varlist [if] [in] , censored(varname)
[level(passthru) *]
 10.                 marksample touse
 11.                 preserve
 12.                 qui keep if `touse'
 13.                 gettoken depvar indepvars : varlist
 14. *Rescale the independent variable into categories one through maxvalue
.                 tempvar newdepvar
 15.                 egen `newdepvar'=group(`depvar')
 16.                 qui sum `newdepvar'
 17.                 local maxvalue =r(max)
 18. *Check to make sure there are uncensored obs for each value
.                 forvalues i=1(1)`maxvalue' {
 19.                         tempvar ucsample sumucs
 20.                         qui gen `ucsample'=0
 21.                         qui replace `ucsample'=(`censored'!=1) if
(`newdepvar'==`i')
 22.                         qui gen `sumucs'=sum(`ucsample')
 23.                         local errorcount =0
 24.                         capture assert `sumucs'[_N]>0&`sumucs'[_N]<.
 25.                         if _rc==9 {
 26.                                 di as error "WARNING: value
number `i' has no uncensored obs"
 27.                         }
 28.                         drop `sumucs'
 29.                         drop `ucsample'
 30.                 }
 31. *Store in local k the number of cutpoints that will be necessary
.                 local k =`maxvalue'-1
 32. *Pass the name of the censoring variable to the evaluator (as a
global macro)
.                 global OGLFXH "`censored'"
 33.                 /*Passes the name of the censoring variable*/
. *Estimate a regular ordered probit to get initial value estimates
.                 qui oprobit `newdepvar' `indepvars'
 34. *Take the coefficients and cutpoints from the ordered probit and
store them.
.                 tempname starting
 35.                 matrix `starting'=e(b)
 36.                 matrix list `starting'
 37.                 tempname cutvector
 38.                 matrix `cutvector'=`starting'[1,"_cut1".."_cut`k'"]
 39. *Put together local macros with the starting values
.                 forvalues i=1(1)`k' {
 40.                         local cut`i'=`cutvector'[1,`i']
 41.                         local initlist "`initlist' /_cut`i'=`cut`i''"
 42.                         local cutlist "`cutlist' /_cut`i'"
 43.                 }
 44.                 foreach var of local indepvars {
 45.                         local `var'_b
=`starting'[1,colnumb(`starting',"`var'")]
 46.                 }
 47.                 foreach var of local indepvars {
 48.                         local initlist "`initlist' `var'=``var'_b'"
 49.                 }
 50. *And finally, estimate the model.
.                 ml model lf coprbt_ll (`newdepvar'=`indepvars' ,
nocons ) `cutlist' ///
>                         if `touse' , maximize title("Homemade Censored Ordered Probit") ///
>                         init(`initlist') `options'
 51.                 ereturn local cmd "coprbt"
 52.         }
 53.         ml display , `level'
 54. end

. *END coprbt.ado
.
. *START coprbt_ll.ado
. *! v. 1.0 for use with coprbt.ado  -Owen Haaga 22 August 2003
. capture program drop coprbt_ll

. program define coprbt_ll
  1.         local i=0
  2.         while "``i''"!="" {
  3.                 local ++i
  4.         }
  5.         local --i                               /*i should
contain the number of arguments*/
  6.         local maxvalue=`i'-1    /*2 other arguments plus k
cutpoints is i, i-1=k+1=maxvalue*/
  7.         local k=`maxvalue'-1    /*k is the total number of cutpoints*/
  8.         forvalues i=1(1)`k' {
  9.                 local cuts "`cuts' cut`i'"
 10.         }
 11.         args lnf theta `cuts'
 12.         tempname cutvector
 13.         matrix `cutvector' =J(1,`k',0)
 14.         forvalues i=1(1)`k' {
 15.                 local cutvalue=`cut`i''[1]
 16.                 /*WEAK.  What if the first obs is not in the sample?*/
.                 matrix `cutvector'[1,`i']=`cutvalue'
 17.         }
 18.         local censored "$OGLFXH"
 19.         /*Any other way to pass this?*/
.         di "$OGLFXH is in OGLFXH"
 20.         tempvar y1 ymax yother
 21.         qui gen `y1'=($ML_y1==1)
 22.         qui gen `ymax'=($ML_y1==`maxvalue')
 23.         qui gen `yother'=($ML_y1!=1)&($ML_y1!=`maxvalue')
 24.         qui replace `lnf'=ln(norm(`cutvector'[1,$ML_y1]-`theta'))
if `y1'==1&`censored'==0
 25.         qui replace
`lnf'=ln(norm(`cutvector'[1,$ML_y1]-`theta')-norm(`cutvector'[1,$ML_y1-1]-`theta'))if
`yother'==1&`censored'==0
 26.         qui replace
`lnf'=ln(1-norm(`cutvector'[1,$ML_y1-1]-`theta')) if
`ymax'==1&`censored'==0
 27.         qui replace `lnf'=ln(1) if `y1'==1&`censored'==1
 28.         qui replace
`lnf'=ln(1-norm(`cutvector'[1,$ML_y1-1]-`theta')) if
`yother'==1&`censored'==1
 29.         qui replace
`lnf'=ln(1-norm(`cutvector'[1,$ML_y1-1]-`theta')) if
`ymax'==1&`censored'==1
 30.         qui replace `lnf'=ln(norm(`cutvector'[1,$ML_y1]-`theta'))
if `y1'==1&`censored'==-1
 31.         qui replace `lnf'=ln(norm(`cutvector'[1,$ML_y1]-`theta'))
if `yother'==1&`censored'==-1
 32.         qui replace `lnf'=ln(1) if `ymax'==1&`censored'==-1
 33. end

.
. *END coprbt_ll.ado
. drawnorm x1 x2 x3 e,n(10000) clear
  -------------------------------------------------------------------------------------------------------------------
begin drawnorm ---
  - local xeqversion : di "version " string(_caller()) ":"
  - version 8.2
  - local version _caller()
  - gettoken first 0: 0, parse(",")
  - syntax [, n(string) SEED(string) Double CORR(string) COV(string)
CStorage(string) Means(string) SDs(string) CLEAR FORCEPSD TOL(passt
> hru) ]
  - quietly count
  - local curn = r(N)
  - if `"`n'"' != "" {
  = if `"10000"' != "" {
  - confirm integer n `n'
  = confirm integer n 10000
  - if `n' == `curn' {
  = if 10000 == 0 {
    local n
    }
  - }
  - if `"`n'"' == "" {
  = if `"10000"' == "" {
    local nobs = r(N)
    if `nobs' <= 0 {
    error 2000
    }
    if "`clear'" != "" {
    drop _all
    qui set obs `nobs'
    local n `nobs'
    }
    }
  - else {
  - if `n' <= 0 {
  = if 10000 <= 0 {
    error 2000
    }
  - qui count
  - if `n' != r(N) {
  = if 10000 != r(N) {
  - qui des, short
  - if r(changed) & ("`clear'" == "" ) {
  = if r(changed) & ("clear" == "" ) {
    error 4
    }
  - drop _all
  - }
  - local nobs = `n'
  = local nobs = 10000
  - qui set obs `nobs'
  = qui set obs 10000
  - }
  - local 0 "`first'"
  = local 0 "x1 x2 x3 e"
  - syntax newvarlist
  - local k : word count `varlist'
  = local k : word count x1 x2 x3 e
  - if "`seed'" != "" {
  = if "" != "" {
    `xeqversion' set seed `seed'
    }
  - tempname C D L M P S
  - if "`cov'" != "" | "`corr'" != "" {
  = if "" != "" | "" != "" {
    if "`cov'" != "" & "`corr'" != "" {
    dis as err "cov() and corr() " "may not be specified together"
    exit 198
    }
    if `"`corr'"' != "" {
    _m2matrix `C' corr `k' "`corr'" "`cstorage'"
    local Cname `corr'
    }
    else {
    _m2matrix `C' cov `k' "`cov'" "`cstorage'"
    local Cname `cov'
    }
    _checkpd `C', matname(`Cname') check(psd) `forcepsd' `tol'
    if r(npos)==`k' {
    matrix `P' = cholesky(`C')
    }
    else {
    matrix `L' = r(L)
    matrix `D' = r(Ev)
    forvalues i = 1/`k' {
    matrix `D'[1,`i'] = sqrt(max(0,`D'[1,`i']))
    }
    matrix `P' = `L'*diag(`D')
    }
    }
  - else {
  - matrix `P' = I(`k')
  = matrix __000004 = I(4)
  - }
  - if "`means'" != "" {
  = if "" != "" {
    _m2matrix `M' means `k' "`means'"
    }
  - else {
  - matrix `M' = J(1,`k', 0)
  = matrix __000003 = J(1,4, 0)
  - }
  - if `"`sds'"' != "" {
  = if `""' != "" {
    if `"`cov'"' != "" {
    dis as err "cov() and sds() may not be specified together"
    exit 198
    }
    _m2matrix `S' sds `k' "`sds'"
    }
  - else {
  - matrix `S' = J(1,`k',1)
  = matrix __000005 = J(1,4,1)
  - }
  - tokenize `varlist'
  = tokenize x1 x2 x3 e
  - local newlist `varlist'
  = local newlist x1 x2 x3 e
  - foreach var of local newlist {
  - if `version' <= 10 {
  = if _caller() <= 10 {
    qui gen `double' `var' = invnormal(uniform())
    }
  - else {
  - qui gen `double' `var' = rnormal()
  = qui gen  x1 = rnormal()
  - }
  - }
  - if `version' <= 10 {
  = if _caller() <= 10 {
    qui gen `double' `var' = invnormal(uniform())
    }
  - else {
  - qui gen `double' `var' = rnormal()
  = qui gen  x2 = rnormal()
  - }
  - }
  - if `version' <= 10 {
  = if _caller() <= 10 {
    qui gen `double' `var' = invnormal(uniform())
    }
  - else {
  - qui gen `double' `var' = rnormal()
  = qui gen  x3 = rnormal()
  - }
  - }
  - if `version' <= 10 {
  = if _caller() <= 10 {
    qui gen `double' `var' = invnormal(uniform())
    }
  - else {
  - qui gen `double' `var' = rnormal()
  = qui gen  e = rnormal()
  - }
  - }
  - mat roweq `P' = " "
  = mat roweq __000004 = " "
  - mat coleq `P' = " "
  = mat coleq __000004 = " "
  - mat rownames `P' = `varlist'
  = mat rownames __000004 = x1 x2 x3 e
  - mat colnames `P' = `varlist'
  = mat colnames __000004 = x1 x2 x3 e
  - forvalues i = 1 / `k' {
  = forvalues i = 1 / 4 {
  - tempname new`i' row
  = tempname new1 row
  - mat `row' = `P'[`i', 1...]
  = mat __000007 = __000004[1, 1...]
  - mat score `new`i'' = `row'
  = mat score __000006 = __000007
  - }
  - tempname new`i' row
  = tempname new2 row
  - mat `row' = `P'[`i', 1...]
  = mat __000009 = __000004[2, 1...]
  - mat score `new`i'' = `row'
  = mat score __000008 = __000009
  - }
  - tempname new`i' row
  = tempname new3 row
  - mat `row' = `P'[`i', 1...]
  = mat __00000B = __000004[3, 1...]
  - mat score `new`i'' = `row'
  = mat score __00000A = __00000B
  - }
  - tempname new`i' row
  = tempname new4 row
  - mat `row' = `P'[`i', 1...]
  = mat __00000D = __000004[4, 1...]
  - mat score `new`i'' = `row'
  = mat score __00000C = __00000D
  - }
  - tokenize `varlist'
  = tokenize x1 x2 x3 e
  - forvalues i = 1 / `k' {
  = forvalues i = 1 / 4 {
  - qui replace ``i'' = `new`i'' * `S'[1,`i'] + `M'[1,`i']
  = qui replace x1 = __000006 * __000005[1,1] + __000003[1,1]
  - }
  - qui replace ``i'' = `new`i'' * `S'[1,`i'] + `M'[1,`i']
  = qui replace x2 = __000008 * __000005[1,2] + __000003[1,2]
  - }
  - qui replace ``i'' = `new`i'' * `S'[1,`i'] + `M'[1,`i']
  = qui replace x3 = __00000A * __000005[1,3] + __000003[1,3]
  - }
  - qui replace ``i'' = `new`i'' * `S'[1,`i'] + `M'[1,`i']
  = qui replace e = __00000C * __000005[1,4] + __000003[1,4]
  - }
  - if "`n'" != "" {
  = if "10000" != "" {
  - dis as txt "(obs `nobs')"
  = dis as txt "(obs 10000)"
(obs 10000)
  - }
  ---------------------------------------------------------------------------------------------------------------------
end drawnorm ---

. g xb=(x1+x2+x3)/2

. g y=(xb+e<-2) +2*inrange(xb+e,-2,0)+3*(xb+e>0)

. tabulate y

          y |      Freq.     Percent        Cum.
------------+-----------------------------------
          1 |        674        6.74        6.74
          2 |      4,380       43.80       50.54
          3 |      4,946       49.46      100.00
------------+-----------------------------------
      Total |     10,000      100.00

. g c=uniform()>.9

. replace y=1 if c==1
(947 real changes made)

.
. coprbt y x1 x2 x3, censored(c)
  ---------------------------------------------------------------------------------------------------------------------
begin coprbt ---
  - version 11.0
  - if replay() {
    if "`e(cmd)'" != "coprbt" {
    error 301
    }
    syntax [, level(passthru)]
    }
  - else {
  - syntax varlist [if] [in] , censored(varname) [level(passthru) *]
  - marksample touse
  - preserve
  - qui keep if `touse'
  = qui keep if __000000
  - gettoken depvar indepvars : varlist
  - tempvar newdepvar
  - egen `newdepvar'=group(`depvar')
  = egen __000001=group(y)
  - qui sum `newdepvar'
  = qui sum __000001
  - local maxvalue =r(max)
  - forvalues i=1(1)`maxvalue' {
  = forvalues i=1(1)3 {
  - tempvar ucsample sumucs
  - qui gen `ucsample'=0
  = qui gen __000002=0
  - qui replace `ucsample'=(`censored'!=1) if (`newdepvar'==`i')
  = qui replace __000002=(c!=1) if (__000001==1)
  - qui gen `sumucs'=sum(`ucsample')
  = qui gen __000003=sum(__000002)
  - local errorcount =0
  - capture assert `sumucs'[_N]>0&`sumucs'[_N]<.
  = capture assert __000003[_N]>0&__000003[_N]<.
  - if _rc==9 {
    di as error "WARNING: value number `i' has no uncensored obs"
    }
  - drop `sumucs'
  = drop __000003
  - drop `ucsample'
  = drop __000002
  - }
  - tempvar ucsample sumucs
  - qui gen `ucsample'=0
  = qui gen __000004=0
  - qui replace `ucsample'=(`censored'!=1) if (`newdepvar'==`i')
  = qui replace __000004=(c!=1) if (__000001==2)
  - qui gen `sumucs'=sum(`ucsample')
  = qui gen __000005=sum(__000004)
  - local errorcount =0
  - capture assert `sumucs'[_N]>0&`sumucs'[_N]<.
  = capture assert __000005[_N]>0&__000005[_N]<.
  - if _rc==9 {
    di as error "WARNING: value number `i' has no uncensored obs"
    }
  - drop `sumucs'
  = drop __000005
  - drop `ucsample'
  = drop __000004
  - }
  - tempvar ucsample sumucs
  - qui gen `ucsample'=0
  = qui gen __000006=0
  - qui replace `ucsample'=(`censored'!=1) if (`newdepvar'==`i')
  = qui replace __000006=(c!=1) if (__000001==3)
  - qui gen `sumucs'=sum(`ucsample')
  = qui gen __000007=sum(__000006)
  - local errorcount =0
  - capture assert `sumucs'[_N]>0&`sumucs'[_N]<.
  = capture assert __000007[_N]>0&__000007[_N]<.
  - if _rc==9 {
    di as error "WARNING: value number `i' has no uncensored obs"
    }
  - drop `sumucs'
  = drop __000007
  - drop `ucsample'
  = drop __000006
  - }
  - local k =`maxvalue'-1
  = local k =3-1
  - global OGLFXH "`censored'"
  = global OGLFXH "c"
  - qui oprobit `newdepvar' `indepvars'
  = qui oprobit __000001  x1 x2 x3
  - tempname starting
  - matrix `starting'=e(b)
  = matrix __000008=e(b)
  - matrix list `starting'
  = matrix list __000008

__000008[1,5]
      __000001:   __000001:   __000001:       cut1:       cut2:
            x1          x2          x3       _cons       _cons
y1   .33256642   .33002666   .31409146  -1.1333436   .14156848
  - tempname cutvector
  - matrix `cutvector'=`starting'[1,"_cut1".."_cut`k'"]
  = matrix __000009=__000008[1,"_cut1".."_cut2"]
_cut2 not found
    forvalues i=1(1)`k' {
    local cut`i'=`cutvector'[1,`i']
    local initlist "`initlist' /_cut`i'=`cut`i''"
    local cutlist "`cutlist' /_cut`i'"
    }
    foreach var of local indepvars {
    local `var'_b =`starting'[1,colnumb(`starting',"`var'")]
    }
    foreach var of local indepvars {
    local initlist "`initlist' `var'=``var'_b'"
    }
    ml model lf coprbt_ll (`newdepvar'=`indepvars' , nocons )
`cutlist' if `touse' , maximize title("Homemade Censored Ordered
Probit")
> init(`initlist') `options'
    ereturn local cmd "coprbt"
    }
  -----------------------------------------------------------------------------------------------------------------------
end coprbt ---
r(111);

end of do-file

r(111);


 I hope now I am able to convey a better picture of what I am doing.
Any help with what is generating the error would be helpful.

Austin told me that for uncensored observations the ordered probit
outperforms the censored ordered probit results. I would still need to
estimate the censored ordered probit results for my data and compare
with what Austin suggested.
Thanks so much.

Urmi

On Mon, Sep 26, 2011 at 1:10 PM, Urmi Bhattacharya <ub3@indiana.edu> wrote:
> Hi Nick,
>
> Thanks. I did not put the input code up because it was taking up too
> much space. I used a simulated data set. I am currently away frrom
> Stata, but I'll put the input code and the simulated data in case that
> helps.
>
> Thank you so much.
>
> Urmi
>
> On Mon, Sep 26, 2011 at 1:02 PM, Nick Cox <njcoxstata@gmail.com> wrote:
>> No. The rest of what he said.
>>
>> But even apart from that, we can't see your input code and dataset.
>> Remote debugging is going to be very very difficult if you don't work
>> with code and data (e.g. simulated data, some set you can -webuse- or
>> -sysuse-) that we can see, not that it will be easy otherwise.
>>
>> Nick
>>
>> On Mon, Sep 26, 2011 at 8:25 AM, Urmi Bhattacharya <ub3@indiana.edu> wrote:
>>> Hi Nick,
>>>
>>> Did you mean the one about lines not being properly structured?
>>>
>>> Thanks
>>>
>>> Urmi
>>>
>>>
>>> On Mon, Sep 26, 2011 at 12:44 PM, Nick Cox <njcoxstata@gmail.com> wrote:
>>>> As far as I can tell you are ignoring Austin's advice in his last email.
>>>>
>>>> I can't illuminate what is happening with your code and data.
>>>>
>>>> Nick
>>>>
>>>> On Mon, Sep 26, 2011 at 3:50 AM, Urmi Bhattacharya <ub3@indiana.edu> wrote:
>>>>> Hi Nick and Austin,
>>>>>
>>>>> Thank you so much for your suggestions. I put the
>>>>> set trace on
>>>>> set traced 1
>>>>>
>>>>> followed by Owen's code and the simulation. It identified some errors
>>>>> I could correct by áwrapping up the unwrapped lines . But I am getting
>>>>> an error now and I can't understand what is causing it. I am providing
>>>>> the result that I am getting
>>>>>
>>>>> - version 11.0
>>>>> á- if replay() {
>>>>> á áif "`e(cmd)'" != "coprbt" {
>>>>> á áerror 301
>>>>> á á}
>>>>> á ásyntax [, level(passthru)]
>>>>> á á}
>>>>> á- else {
>>>>> á- syntax varlist [if] [in] , censored(varname) [level(passthru) *]
>>>>> á- marksample touse
>>>>> á- preserve
>>>>> á- qui keep if `touse'
>>>>> á= qui keep if __000000
>>>>> á- gettoken depvar indepvars : varlist
>>>>> á- tempvar newdepvar
>>>>> á- egen `newdepvar'=group(`depvar')
>>>>> á= egen __000001=group(y)
>>>>> á- qui sum `newdepvar'
>>>>> á= qui sum __000001
>>>>> á- local maxvalue =r(max)
>>>>> á- forvalues i=1(1)`maxvalue' {
>>>>> á= forvalues i=1(1)3 {
>>>>> á- tempvar ucsample sumucs
>>>>> á- qui gen `ucsample'=0
>>>>> á= qui gen __000002=0
>>>>> á- qui replace `ucsample'=(`censored'!=1) if (`newdepvar'==`i')
>>>>> á= qui replace __000002=(c!=1) if (__000001==1)
>>>>> á- qui gen `sumucs'=sum(`ucsample')
>>>>> á= qui gen __000003=sum(__000002)
>>>>> á- local errorcount =0
>>>>> á- capture assert `sumucs'[_N]>0&`sumucs'[_N]<.
>>>>> á= capture assert __000003[_N]>0&__000003[_N]<.
>>>>> á- if _rc==9 {
>>>>> á ádi as error "WARNING: value number `i' has no uncensored obs"
>>>>> á á}
>>>>> á- drop `sumucs'
>>>>> á= drop __000003
>>>>> á- drop `ucsample'
>>>>> á= drop __000002
>>>>> á- }
>>>>> á- tempvar ucsample sumucs
>>>>> á- qui gen `ucsample'=0
>>>>> á= qui gen __000004=0
>>>>> á- qui replace `ucsample'=(`censored'!=1) if (`newdepvar'==`i')
>>>>> á= qui replace __000004=(c!=1) if (__000001==2)
>>>>> á- qui gen `sumucs'=sum(`ucsample')
>>>>> á= qui gen __000005=sum(__000004)
>>>>> á- local errorcount =0
>>>>> á- capture assert `sumucs'[_N]>0&`sumucs'[_N]<.
>>>>> á= capture assert __000005[_N]>0&__000005[_N]<.
>>>>> á- if _rc==9 {
>>>>> á ádi as error "WARNING: value number `i' has no uncensored obs"
>>>>> á á}
>>>>> á- drop `sumucs'
>>>>> á= drop __000005
>>>>> á- drop `ucsample'
>>>>> á= drop __000004
>>>>> á- }
>>>>> á- tempvar ucsample sumucs
>>>>> á- qui gen `ucsample'=0
>>>>> á= qui gen __000006=0
>>>>> á- qui replace `ucsample'=(`censored'!=1) if (`newdepvar'==`i')
>>>>> á= qui replace __000006=(c!=1) if (__000001==3)
>>>>> á- qui gen `sumucs'=sum(`ucsample')
>>>>> á= qui gen __000007=sum(__000006)
>>>>> á- local errorcount =0
>>>>> á- capture assert `sumucs'[_N]>0&`sumucs'[_N]<.
>>>>> á= capture assert __000007[_N]>0&__000007[_N]<.
>>>>> á- if _rc==9 {
>>>>> á ádi as error "WARNING: value number `i' has no uncensored obs"
>>>>> á á}
>>>>> á- drop `sumucs'
>>>>> á= drop __000007
>>>>> á- drop `ucsample'
>>>>> á= drop __000006
>>>>> á- }
>>>>> á- local k =`maxvalue'-1
>>>>> á= local k =3-1
>>>>> á- global OGLFXH "`censored'"
>>>>> á= global OGLFXH "c"
>>>>> á- qui oprobit `newdepvar' `indepvars'
>>>>> á= qui oprobit __000001 áx1 x2 x3
>>>>> á- tempname starting
>>>>> á- matrix `starting'=e(b)
>>>>> á= matrix __000008=e(b)
>>>>> á- tempname cutvector
>>>>> á- matrix `cutvector'=`starting'[1,"_cut1".."_cut`k'"]
>>>>> á= matrix __000009=__000008[1,"_cut1".."_cut2"]
>>>>>
>>>>> _cut2 not found
>>>>>
>>>>> á áforvalues i=1(1)`k' {
>>>>> á álocal cut`i'=`cutvector'[1,`i']
>>>>> á álocal initlist "`initlist' /_cut`i'=`cut`i''"
>>>>> á álocal cutlist "`cutlist' /_cut`i'"
>>>>> á á}
>>>>> á áforeach var of local indepvars {
>>>>> á álocal `var'_b =`starting'[1,colnumb(`starting',"`var'")]
>>>>> á á}
>>>>> á áforeach var of local indepvars {
>>>>> á álocal initlist "`initlist' `var'=``var'_b'"
>>>>> á á}
>>>>> á áml model lf coprbt_ll (`newdepvar'=`indepvars' , nocons ) `cutlist'
>>>>> á áif `touse' , maximize title("Homemade Censored Ordered Probit")
>>>>> init(`initlist') `options'
>>>>> á áereturn local cmd "coprbt"
>>>>> á á}
>>>>>
>>>>> r(111);
>>>>>
>>>>> end of do-file
>>>>>
>>>>> r(111);
>>>>>
>>>>> The "_cut2 not found" is in red suggesting that is where the problem
>>>>> is, but I cannot figure out what is causing it.
>>>>>
>>>>> Any suggestions?
>>>>>
>>>>> Thanks a lot.
>>>>>
>>>>> Urmi
>>>>> On Sat, Sep 24, 2011 at 8:14 PM, Austin Nichols <austinnichols@gmail.com> wrote:
>>>>>> Urmi Bhattacharya <ub3@indiana.edu>:
>>>>>> I also guess that you have not got the line structure right in your
>>>>>> files, but I will also note that a quick simulation I ran indicated
>>>>>> that -oprobit- run only on uncensored obs outperforms unrestricted
>>>>>> -oprobit- and -coprbit-. áSo start with
>>>>>>
>>>>>> . áoprobit y x1 x2 x3 if censored==0
>>>>>>
>>>>>> I am also guessing that other models would outperform all 3 options
>>>>>> mentioned here, but that will have to wait for another day.
>>>>>>
>>>>
>>
>> *
>> *   For searches and help try:
>> *   http://www.stata.com/help.cgi?search
>> *   http://www.stata.com/support/statalist/faq
>> *   http://www.ats.ucla.edu/stat/stata/
>>
>

*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/statalist/faq
*   http://www.ats.ucla.edu/stat/stata/


© Copyright 1996–2014 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   Site index