Probit Analysis 1


Overview

Version info: Code for this page was tested in Gretl 1.9.92

 Probit regression, also called a probit model, is used to model dichotomous or binary outcome variables. In the probit model, the inverse standard normal distribution of the probability is modeled as a linear combination of the predictors. 

 Please Note: The purpose of this page is to show how to use various data analysis commands. It does not cover all aspects of the research process which researchers are expected to do. In particular, it does not cover data cleaning and checking, verification of assumptions, model diagnostics and potential followup analyses.

Description of the data

For our data analysis below, we are going to use an example about getting into graduate school. Hypothetical data was generated by UCLA, which can be obtained from this website.
http://www.ats.ucla.edu/stat/stata/dae/binary.dta

 However, in gretl you can use the command log to access it, and that's the next step: 
? open http://www.ats.ucla.edu/stat/data/binary.csv

parsing C:\Users\Jaysp_000\AppData\Roaming\gretl\binary.csv...

using delimiter ','

     longest line: 20 characters 

     first field: 'admit'

     number of columns = 4

     number of variables: 4

     number of non-blank lines: 401

scanning for variable names...

     line: admit,gre,gpa,rank
scanning for row labels and data...

treating these as undated data



Listing 5 variables:
  0) const    1) admit    2) gre      3) gpa      4) rank

This data set has a binary response (outcome, dependent) variable called admit. There are three predictor variables: gre, gpa and rank. We will treat the variables gre and gpa as continuous. The variable rank is ordinal, it takes on the values 1 through 4. Institutions with a rank of 1 have the highest prestige, while those with a rank of 4 have the lowest. We will treat rank as categorical.

 The odd thing with Gretl is that even though it can recognize categorical variables, it's not going to automatically convert it for you. Right now, we are dealing with a single variable in rank that so happens to carry the mix of discrete values 1 to 4. So in order to turn rank into a categorical variable in Gretl, we need to first convert it to (or just make sure it is) a discrete variable. In order to do that, we dummify it. No really, dummify it. Dummify (or 'dummify' more precisely) is a function in Gretl used to turn discrete variables into categorical ones. For every categorical variable in a regression estimation, there is always a baseline, or reference group, to calculate partial effects. I will get into that later, but for now...

? setinfo rank --discrete
? dummify rank --drop-first #reference group
Listing 8 variables:
     0) const       1) admit    2) gre       3) gpa     4) rank
     5) Drank_2     6) Drank_3      7) Drank_4


Now we have a list of 7 variables to work with: the original set, and the extra variables. Notice there are only 3 extra variables (for ranks 2 to 4) and rank 1 is not included. This is because rank 1 is the reference group, or the base that we analyze the partial effects from. Again, we will see more this, but for now it should be mentioned that rank could have easily been split into all 4 category variables, but I decided against it by skipping the first rank using the '--drop-first' option along with the dummify command. In another scenario, if you wanted to drop the last variable, you can just use the '--drop-last' option instead. 

Before we continue, we should probably do something about the three new dummy variables. The names should be a little more convenient, just in case we need it. Notice how each variable (particularly Drank_2 to Drank_4) have corresponding numbers reflecting their position on the list of variables. This is necessary as we use the 'rename' command.
? rename 5 rank2  #replacing the name at 5th spot with 'rank2'
? rename 6 rank3  #replacing the name at 6th spot with 'rank3'
? rename 7 rank4  #replacing the name at 7th spot with 'rank4'

Moving forward: Summary statistics and more

It is very much good practice to check with the basic statistics before you regress the variables.

#Generate a list so I can summarize all the data 
? list mydata = admit gre gpa rank
Generated list mydata

#Now summarize the data
? summary mydata

                     Mean         Median        Minimum        Maximum
admit             0.31750        0.00000        0.00000         1.0000
gre                587.70         580.00         220.00         800.00
gpa                3.3899         3.3950         2.2600         4.0000
rank               2.4850         2.0000         1.0000         4.0000

                Std. Dev.           C.V.       Skewness   Ex. kurtosis
admit             0.46609         1.4680        0.78410        -1.3852
gre                115.52        0.19656       -0.14381       -0.34202
gpa               0.38057        0.11226       -0.20791       -0.58713
rank              0.94446        0.38006       0.096858       -0.89544

                 5% perc.      95% perc.       IQ range   Missing obs.
admit             0.00000         1.0000         1.0000              0
gre                381.00         800.00         140.00              0
gpa                2.7315         4.0000        0.54000              0
rank               1.0000         4.0000         1.0000              0

#Frequency distribution tables for 'rank' and 'admit'
? freq rank

Frequency distribution for rank, obs 1-400

          frequency    rel.     cum.

   1          61     15.25%   15.25% *****
   2         151     37.75%   53.00% *************
   3         121     30.25%   83.25% **********
   4          67     16.75%  100.00% ******

? freq admit

Frequency distribution for admit, obs 1-400

          frequency    rel.     cum.

   0         273     68.25%   68.25% ************************
   1         127     31.75%  100.00% ***********

#Cross-tabulate rank + admit
? xtab rank admit

Cross-tabulation of rank (rows) against admit (columns)

       [   0][   1]  TOT.

[   1]    28    33     61
[   2]    97    54    151
[   3]    93    28    121
[   4]    55    12     67

TOTAL    273   127    400

Pearson chi-square test = 25.2421 (3 df, p-value = 1.37413e-005)

Analysis methods you might consider:

Below is a list of some analysis methods you may have encountered. Some of the methods listed are quite reasonable while others have either fallen out of favor or have limitations.
  • Probit regression, the focus of this page.
  • Logistic regression. A logit model will produce results similar probit regression. The choice of probit versus logit depends largely on individual
  • preferences.
  • OLS regression. When used with a binary response variable, this model is known as a linear probability model and can be used as a way to describe conditional probabilities. However, the errors (i.e., residuals) from the linear probability model violate the homoskedasticity and normality of errors assumptions of OLS regression, resulting in invalid standard errors and hypothesis tests. For a more thorough discussion of these and other problems with the linear probability model, see Long (1997, p. 3840).
  • Two-group discriminant function analysis. A multivariate method for dichotomous outcome variables.

Probit regression

Below the probit command is used to estimate a probit regression model. Because we are working with a categorical variable (aka factor variable) in rank (along with rank2, rank3, rank4), you're going to have to take this series of indicator variables into account. You can just use the rank variable in the regression instead, but it wouldn't help you very much because you don't know if there are any significant partial effects to consider, so it should be included in the model. Note that you can also access this probit estimator from the pull-down menus (on the GUI) using Model>Nonlinear models>Probit>Binary.
? probit admit const gre gpa rank2 rank3 rank4

Model 1: Probit, using observations 1-400
Dependent variable: admit
Standard errors based on Hessian

             coefficient   std. error      z         slope  
  ------------------------------------------------------------
  const      −2.38684      0.674088      −3.541            
  gre         0.00137559   0.000648865    2.120    0.000479833
  gpa         0.477730     0.195463       2.444    0.166642
  rank2      −0.415399     0.195377      −2.126   −0.140263
  rank3      −0.812138     0.208596      −3.893   −0.253574
  rank4      −0.935899     0.245634      −3.810   −0.261630

Mean dependent var   0.317500   S.D. dependent var   0.466087
McFadden R-squared   0.083131   Adjusted R-squared   0.059129
Log-likelihood      −229.2066   Akaike criterion     470.4132
Schwarz criterion    494.3620   Hannan-Quinn         479.8972

Number of cases 'correctly predicted' = 284 (71.0%)
f(beta'x) at mean of independent vars = 0.349
Likelihood ratio test: Chi-square(5) = 41.5633 [0.0000]

           Predicted
              0     1
  Actual 0  254    19
         1   97    30

Test for normality of residual -
  Null hypothesis: error is normally distributed
  Test statistic: Chi-square(2) = 0.286696
  with p-value = 0.866452

  • At the top of the output we see that all 400 observations in our data set were used in the analysis (fewer observations would have been used if any of our variables had missing values).
  • The likelihood ratio chi-square of 41.56 with a p-value of 0.0000 in gretl, but 0.0001 on Stata. It's just to show that gretl often times wont round the number off as opposed to other software. Either way, it tells us that our model as a whole is statistically significant, that is, it fits significantly better than a model with no predictors. 
  • In the table we see the coefficients, their standard errors, the z-statistic, and associated slope (which captures the marginal effects at the means) along with the measures of fit table for statistics.
  • Automatically printed at the bottom is a test showing us that the residuals are normally distributed and just above that are the results of the how much cases in the model was predicted right vs. wrong.
    •
If you're still on the command line, and you wanted to see the p-values, one of the methods you could just use the '--p-values' option right after the regression equation.

? probit admit const gre gpa rank2 rank3 rank4 --p-values

Model 2: Probit, using observations 1-400
Dependent variable: admit
Standard errors based on Hessian

             coefficient   std. error      z      p-value
  --------------------------------------------------------
  const      −2.38684      0.674088      −3.541   0.0004   ***
  gre         0.00137559   0.000648865    2.120   0.0340   **
  gpa         0.477730     0.195463       2.444   0.0145   **
  rank2      −0.415399     0.195377      −2.126   0.0335   **
  rank3      −0.812138     0.208596      −3.893   9.89e-05 ***
  rank4      −0.935899     0.245634      −3.810   0.0001   ***

For convenience, I cut off the rest of the results: they were already showed in the previous model.
  • Both gre, gpa, and the three indicator variables for rank are statistically significant. The probit regression coefficients give the change in the zscore or probit index for a one unit change in the predictor.
  • Notice the p-values replaced with slope (or marginal effects) section.
    • For a one unit increase in gre, the zscore increases by 0.001. 
    • For each one unit increase in gpa, the zscore increases by 0.478. 
    • The indicator variables for rank have a slightly different interpretation. For example, having attended an undergraduate institution of rank of 2, versus an institution with a rank of 1 (the reference group), decreases the zscore by 0.415.
You might notice that the confidence intervals and significance indicators are not present in this current output. If you're estimating using the GUI approach, you can find the confidence intervals by going to Analysis > Confidence intervals for coefficients.

Supplemental tests

We can test for an overall effect of rank using the Wald Test. Below we see that the overall effect of rank is statistically significant. On gretl that is done using the 'omit' command and then '--test-only' command to show the results.

? omit rank2 rank3 rank4 --chi-square --test-only

Test on Model 2:

  Null hypothesis: the regression parameters are zero for the variables
    rank2, rank3, rank4
  Wald test: Chi-square(3) = 21.3169, p-value 9.04679e-005
  (F-form: F(3, 394) = 7.10564, p-value 0.000116592)

The chi-squared test statistic of approximately 21.4 with 3 degrees of freedom is associated with a p-value of less than 0.001, indicating that the overall effect of rank is statistically significant.

Linear Restrictions. We can also test additional hypotheses about the differences in the coefficients for different levels of rank. Below I test that the coefficient for rank=2 is equal to the coefficient for rank=3. 

? restrict
? b[4]-b[5]= 0    # or b[rank2]-b[rank3]=0
? end restrict

Restriction:
 b[rank2] - b[rank3] = 0

Test statistic: chi^2(1) = 5.60168, with p-value = 0.0179433

Predicted Probabilities can also be used to help you understand the model. In Gretl I prefer to calculate this right after the probit estimates (sometimes you can unknowingly alter your data with subsequent calculations). You cannot directly calculate predicted probabilities using any kind of command log or drop down menu in gretl, unlike other software. Below I created a simple function to calculate a simple predicted probability. The parameters needed are the hypothetical gre, gpa, rank (from the rank2 to rank4 dummies), and the coefficients, all in order. Lets take on three examples:

  1. The predicted probability of admission if a person went to a school that is ranked as least prestigious (lets use 4), has a GPA of 3.7, and has a GRE score of 1800.
  2. The predicted probability of admission if a person went to the highest ranked school (lets say rank 1, has a 3.0 GPA and has a GRE score of 1440. 
  3. The predicted probability of admission if a person went to a relatively high ranked school (lets use 2), has a 4.0 GPA, and has a GRE score of 1300.
# Just to be safe, we regress again, but keep the results hidden this time
? probit admit const gre gpa rank2 rank 3 rank 4 --quiet

# create the PP function
? function scalar predictor (scalar GRE, scalar GPA, scalar rank2, scalar rank3, scalar rank4, matrix param)
> matrix xx ={1, GRE, GPA, rank2, rank3, rank4}
> scalar p = xx*param
> return cnorm(p)
> end function

## Example 1:
? matrix ex1 = predictor (1800,3.7,0,0,1,$coeff)
? ex1
Generated scalar ex1 = 0.82146

## Example 2:
? matrix ex2 = predictor (1440,3,0,0,0,$coeff)
? print ex2
Generated scalar ex2 = 0.84784

## Example 3:
? matrix ex3 = predictor (1300,4,1,0,0,$coeff)
? ex3
Generated scalar ex3 = 0.81513

The 3 examples give us three different predicted probabilities. A person matching the description in example 1 is about 82% more likely to get admitted, which is slightly trailing the chances of the person in example 2, who goes to a rank 1 school but has a lower GPA and GRE score. That person is approximately 85% more likely of getting admitted, and the person last example is about 82% more likely to be admitted. If anything, we can make a brief inference that there is heavier influence by school prestige in rank than GRE scores when it comes to getting into grad school.

Things to consider

  • Empty cells or small cells: You should check for empty or small cells by doing a cross tabulation between categorical predictors and the outcome variable. If a cell has very few cases (a small cell), the model may become unstable or it might not run at all.
  • Separation or quasi-separation (also called perfect prediction), a condition in which the outcome does not vary at some levels of the independent variables. See our page FAQ: What is complete or quasi-complete separation in logistic/probit regression and how do we deal with them? for information on models with perfect prediction.
  • Sample size: Both probit and logit models require more cases than OLS regression because they use maximum likelihood estimation techniques. It is sometimes possible to estimate models for binary outcomes in datasets with only a small number of cases using exact logistic regression (using the exlogistic command). For more information see our data analysis example for exact logistic regression. It is also important to keep in mind that when the outcome is rare, even if the overall dataset is large, it can be difficult to estimate a probit model.
  • PseudoRsquared: Many different measures of psuedo-Rsquared exist. One of these measures, the McFadden R-squared, is provided in Gretl, in the measures of fit section when you run a regression estimation. They all attempt to provide information similar to that provided by Rsquared in OLS regression; however, none of them can be interpreted exactly as Rsquared in OLS regression is interpreted. For a discussion of various pseudo-Rsquareds see Long and Freese (2006) or UCLAs FAQ page What are pseudo Rsquareds? 
  • Diagnostics: The diagnostics for probit regression are different from those for OLS regression. The diagnostics for probit models are similar to those for logit models. For a discussion of model diagnostics for logistic regression, see Hosmer and Lemeshow (2000, Chapter 5).

Sources:

Heavily borrowed and influenced from:
  • Gretl User's Guide: Jan 2015
  • Hosmer, D. & Lemeshow, S. (2000). Applied Logistic Regression (Second Edition). New York: John Wiley & Sons, Inc.
  • Long, J. Scott (1997). Regression Models for Categorical and Limited Dependent Variables. Thousand Oaks, CA: Sage Publications.
  • Tutorial page from UCLA

Joshua

The logic of Binary Dependent Variable Modeling

In an regression, you're often interested in effects. Just how does the variable of interest (Y) get affected by the change in variables which are assumed to be independent (Xs). 

Regression with continuous dependent variables

Lets say you're all of a sudden interested in the attendance to baseball games. You want to understand what affects the level of attendance to baseball games, lets call it ATTEND (in thousands of people), and you believe that a current wins in the season have a say in how the level of attendance to the baseball games look. In other words, you believe there is a linear relationship between attendance to a baseball team's games and the wins they have in a season. In that case, we can notate the model like this:


 or simply...

 Based off this linear relationship, estimating this model is fairly easy. You can observe the relationship between the two using graphical charts like scatter plots, and based on the apparent upward trend in the plot points, you may infer that there is a positive relationship between the X and Y variables, probably leaving you with the question of how strong and significant the magnitude of this relationship is.

With this particular example, you're bound to receive regression estimates that confirm the positive relationship between the two variables.

Dependent variable: ATTEND

              coefficient    std. error   t-ratio   p-value
  -----------------------------------------------------------
  const        −1340.92      387.127      −3.464    0.0009    ***
  CURNTWIN      38.6325      4.74786       8.137    6.02e-012 ***

And when you interpret the magnitude of the coefficients, you can conclude that for every unit increase in CURNTWIN, ATTEND sees an increase by 38.6325. In laymen's terms, for every additional win a baseball team has, their stadiums are expected to experience an increase in attendance by about approximately 38 to 39 people.

OLS models are very popular approaches to regression analysis, and in these models we're dealing with dependent variables that have continuous values.

Regression with Dichotomous Dependent variables

Lets change that up and assume that instead of baseball wins and attendances, you want to find out the probability that a person attends college given the percentage level of their parental wage. We can just write that as...
 or simply...

We're going to move away from the linear combination of the independent variables (LCIV) because they actually give off nonsensical results. For instance, a regression equation with a continuous variable like income can provide useful information when you generate graphs, like a scatter plot, but with dichotomous variable as Y, we get something where all the plot points are scattered to either one value or the another, like this:

The points are bounded to either 0 or 1, which makes a lot of sense (if you're familiar with dummy variables) because 0 and 1 represent outcomes (0 ="didn't attend college" ; 1 = "attended college").

So the rectification is to take a non-linear transformation of the linear combination of independent variables, looking like this:


The properties we expect function F to have reflect the discrete outcomes of interest (as related to qualitative variables). In that role, function F is expected to output 0 as the linear combination of independent variables tend towards negative infinity, and F is expected to output 1 as the LCIVs tend towards positive infinity. That's just a mathematical way of stating this hypothetical point, that as we see a continual increase (or decrease) in the value(s) we get from the LCIVs, we should also gradually expect the outcome of interest (Y) to be a likely occurrence (or likely nonoccurrence). 

The corresponding math notations usually looks like this:

F(-∞) = 0
F(+∞) = 1

and thus, this transformation of the independent variables means that our probability belongs in the range of 0 and 1, or...

In context, lets assume that we see the value for X as very high, indicating high parental income; in that case, then we can assume that the probability of a person going to college is virtually 1, or very likely to happen.

Conclusions

The basic logic of logit and probit modeling is grounded on this reasoning. Instead of trying to see the effects of the X var(s) on the continuous Y variable, like the effects of income on overall consumption, you're trying to see the how the X var(s) influences the probability that dichotomous Y variable will occur; in other words, the likelihood of outcomes.

 It should be understood and reminded that yes, we can expect a higher or lower probability of an event occurring but that does not mean it will indeed happen or not. The ultimate outcome can most definitely go the other way. 

 It should also be noted that regression equations always include error terms, which were purposely excluded in this post. Just know that you are to always include the error term when conducting analyses for academic or professional work.

Sources:

Heavily borrowed and influenced from:

Joshua

How to Analyze a Business Model

A post by digital marketing professional, Eric Noren, guides the reader on how to analyze business models. I find this post very useful and informative. Its not too long, and summarizes a bunch of the basics. Click here (or the picture) to read it.

Joshua

Subscribe to: Posts (Atom)