Simple linear and nonlinear models using Stata’s ml command | Stata Code Fragments

This code fragment page is designed to show how a number of simple linear and nonlinear models can be programmed using Stata’s ml command. These programs are not complete estimation commands but just those parts needed to compute the relevent statistics.

Acknowlegements: Most of the code on this page was borrowed or adapted from various sources that I would like to acknowledge here.

The code for OLS, binary logistic and probit regression came from “Maximum Likelihood Estimation with Stata,” by William Gould, Jeffrey Pitblado, and William Sribney

The code for poisson and negative binomial regression came from “Microeconometrics Using Stata,” by A. Colin Cameron and Pravin K. Trivedi.

The code for multinomial logistic regression was adapted from a Statalis posting by Maarten Buis.

All the programs use the same dataset, hsbdemo, and the same predictor variables, read and female. Please note that the programs for ordered logistic, ordered probit and multinomial logistic regression are coded to work with variables having exactly three categories. If your variables have more than three categories the code will need to be modified.

One way to tryout these programs is to copy and paste all of the code into your do-file editor. Load the hsbdemo dataset. Then, highlight the sections of code you are interested in and run it by clicking on the “Do document/selection” button in the editor.

Here are the code fragments:.

use https://stats.idre.ucla.edu/stat/data/hsbdemo


/****************************/
/*    linear regression     */
/****************************/

capture program drop lfols
program lfols
  version 10.1
  args lnf xb lnsigma
  local y "$ML_y1"
  quietly replace `lnf' = ln(normalden(`y', `xb',exp(`lnsigma')))
end 

 ml model lf lfols (xb: write = read female) (lnsigma:)
 ml maximize
 display exp([lnsigma]_cons)

/**********************************/
/*   binary logistic regression   */
/**********************************/

capture program drop lflogit
program lflogit
  version 10.1
  args lnf xb
  local y "$ML_y1"
  quietly replace `lnf' = ln(  invlogit(`xb')) if `y'==1
  quietly replace `lnf' = ln(1-invlogit(`xb')) if `y'==0
end

 ml model lf lflogit (honors = read female)
 ml maximize


/*************************/
/*   probit regression   */
/*************************/

capture program drop lfprobit
program lfprobit
  version 10.1
  args lnf xb
  local y "$ML_y1"
  quietly replace `lnf' = ln(  normal(`xb')) if `y'==1
  quietly replace `lnf' = ln(1-normal(`xb')) if `y'==0
end

 ml model lf lfprobit (honors = read female)
 ml maximize

/*************************************/
/*   ordered logitistic regression   */
/*************************************/

capture program drop lfologit
program lfologit
  version 10.1
  args lnf xb a1 a2
  local y "$ML_y1"
  quietly replace `lnf' =  ln(  invlogit(`a1'-`xb')) if `y'==1
  quietly replace `lnf' =  ln(  invlogit(`a2'-`xb') -invlogit(`a1'-`xb'))  if `y'==2
  quietly replace `lnf' =  ln(1-invlogit(`a2'-`xb')) if `y'==3
end

 ml model lf lfologit (ses = read female, nocons) /cut1 /cut2
 ml maximize
 
/*********************************/
/*   ordered probit regression   */
/*********************************/

capture program drop lfoprobit
program lfoprobit
  version 10.1
  args lnf xb a1 a2
  local y "$ML_y1"
  quietly replace `lnf' =  ln(  normal(`a1'-`xb')) if `y'==1
  quietly replace `lnf' =  ln(  normal(`a2'-`xb') - normal(`a1'-`xb'))  if `y'==2
  quietly replace `lnf' =  ln(1-normal(`a2'-`xb')) if `y'==3
end

 ml model lf lfoprobit (ses = read female, nocons) /cut1 /cut2
 ml maximize
 
/***************************************/
/*   multinomial logistic regression   */
/***************************************/


tab prog, gen(pr)
global y1 "pr1"
global y2 "pr2"
global y3 "pr3"

capture program drop lfmlogit
program lfmlogit
  version 10.1
  args lnf xb1 xb2
  tempvar p1 p2 p3
  quietly {
     gen double `p1' = 1/(1+exp(`xb1')+exp(`xb2'))
     gen double `p2' = exp(`xb1')/(1+exp(`xb1')+exp(`xb2'))
     gen double `p3' = exp(`xb2')/(1+exp(`xb1')+exp(`xb2'))
     replace `lnf' = $y1*ln(`p1') + $y2*ln(`p2') + $y3*ln(`p3')
   }
end

 ml model lf lfmlogit (eq1: prog = read female) (eq2: read female) 
 ml maximize

/**************************/
/*   poisson regression   */
/**************************/

capture program drop lfpois
program lfpois
  version 10.1
  args lnf theta1
  tempvar lnyfact mu
  local y "$ML_y1"
  generate double `lnyfact' =lnfactorial(`y')
  generate double `mu' = exp(`theta1')
  quietly replace `lnf' = -`mu' + `y'*`theta1' - `lnyfact'
 end
 
 ml model lf lfpois (awards = read female), vce(robust)
 ml maximize
 
/************************************/ 
/*   negative binomial regression   */
/************************************/
 
 capture program drop lfnb
 program lfnb
   version 10.1
   args lnf theta1 a
   tempvar mu
   local y "$ML_y1"
   generate double `mu' = exp(`theta1')
   quietly replace `lnf' = lngamma(`y'+ (1/`a'))-lngamma((1/`a')) ///
           - lnfactorial(`y') - (`y'+(1/`a'))*ln(1+`a'*`mu')      ///
           + `y'*ln(`a') + `y'*ln(`mu')
 end
 
 ml model lf lfnb (awards = read female)(a:)
 ml maximize