Lecture 27 | Poisson regression

1 STATS 200: introduction to Statistical InferenceAutumn 2016 Lecture 27 Poisson The Poisson log-linear modelExample in the central nervous system transmit signals via a series of actionpotentials, or spikes . The spiking of a single neuron may be measured by a microelectrode,and its sequence of spikes over time is called a spike train. A simple and commonly-usedstatistical model for a spike train is an inhomogeneous Poisson point process, which hasthe following property: Forntime windows of length , lettingYidenote the number ofspikes generated by the neuron in theithtime window, the random variablesY1,..,Ynareindependent and distributed asYi Poisson ( i ), where the parameter icontrols thespiking rate in theithtime window.

2 For simplicity, we will assume = spiking rate iof a neuron may be influenced by external sensory stimuli present inthisithwindow of time, for example the intensity and pattern of light visible to the eye orthe texture of an object presented to the touch. To understand the effects of these sensorystimuli on the spiking rate of a particular neuron, we may perform an experiment that appliesdifferent stimuli in different windows of time and records the neural response. Encoding thestimuli applied in theithwindow of time by a set ofpcovariatesxi1,..,xip, a simple modelfor the Poisson rate parameter iis given bylog i= 0+ 1xi1+..+ pxip,( )or equivalently, i=e 0+ 1xi1+.

3 + with the distributional assumptionYi Poisson ( i), this is called thePoissonlog-linear model, or the Poisson regression model. It is a special case of what is known inneuroscience as the linear-nonlinear Poisson cascade generally, the Poisson log-linear model is a model fornresponsesY1,..,Ynthattake integer count values. EachYiis modeled as an independent Poisson ( i) random variable,where log iis a linear combination of the covariates corresponding to in the cases of linear and logistic regression , we treat the covariates as fixed constants,and the model parameters to be inferred are the regression coefficients = ( 0,.., p). Statistical inferenceWe will describe the procedure for maximum-likelihood estimation of the regression coeffi-cients and Fisher-information based estimation of their standard errors, and discuss someissues concerning model misspecification and robust standard error.

4 ,Ynare independent Poisson random variables, the likelihood function isgiven bylik( 0,.., p) =n i=1 Yiie iYi!where iis defined in terms of 0,.., pand the covariatesxi1,..,xipvia equation ( ).Settingxi0 1 for alli, the log-likelihood is thenl( 0,.., p) =n i=1 Yilog i i logYi!=n i=1Yi(p j=0 jxij) e pj=0 jxij logYi!and the MLEs are the solutions to the system of score equations, form= 0,..,p,0 = l m=n i=1xim(Yi e pj=0 jxij).These equations may be solved numerically using the Newton-Raphson Fisher information matrixIY( ) = E [ 2l( )] may be obtained by computing thesecond-order partial derivatives ofl: 2l m l= n i=1ximxile pj=0 (x1j.)

5 ,xnj) as thejth column of the covariate matrixXand defining thediagonal matrixW=W( ) := diag(e pj=0 jx1j,..,e pj=0 jxnj),the above may be written as 2l m l= XTmWXl, so 2l( ) = XTWXandIY( ) =XTWX. For largen, if the Poisson log-linear model is correct, then the MLE vector isapproximately distributed asN( ,(XTWX) 1). We may then estimate the standard errorof jby sej= ((XT WX) 1)jj,where W=W( ) is the plugin estimate forW. These formulas are the same as for the caseof logistic regression in Lecture 26, except with a different form of the diagonal modeling assumption of a Poisson distribution forYiis rather restrictive, as it impliesthat the variance ofYimust be equal to its mean.

6 This is rarely true in practice, and it isfrequently the case that the observed variance ofYiis larger than its mean this problem isknown asoverdispersion. Nonetheless, the Poisson regression model is oftentimes used inoverdispersed settings: As long asY1,..,Ynare independent andlogE[Yi] = 0+ 1xi1+..+ pxip27-2for eachi(so the model for the means of theYi s is correct), then it may be shown that theMLE in the Poisson regression model is unbiased for , even if the distribution ofYiisnot Poisson and the variance ofYiexceeds its mean. The above standard error estimate sejand the associated confidence interval for j, though, would not correct in the overdispersedsetting.

7 One may use instead the robust sandwich estimate of the covariance of , given by(XT WX) 1(XT WX)(XT WX) 1where W= diag((Y1 1)2,..,(Yn n)2)and i=e pj=0 jxijis the fitted value of for theithobservation. Alternatively, one mayuse the pairs bootstrap procedure as described in Lecture linear model, logistic regression model, and Poisson regression modelare all examples of thegeneralized linear model (GLM). In a generalized linear model,Y1,..,Ynare modeled as independent observations with distributionsYi f(y| i) for someone-parameter familyf(y| ). The parameter iis modeled asg( i) = 0+ 1xi1+..+ pxipfor some one-to-one transformationg:R Rcalled thelink function, wherexi1.

8 ,xipare covariates corresponding toYi. In the linear model considered in Lecture 25, the pa-rameter was wheref(y| ) was the PDF of theN( , 20) distribution (for a knownvariance 20), andg( ) = . In logistic regression , the parameter was pwheref(y|p)was the PMF of the Bernoulli(p) distribution, andg(p) = logp1 p. In Poisson regression ,the parameter was wheref(y| ) was the PMF of the Poisson ( ) distribution, andg( ) = log .The choice of the link functiongis an important modeling decision, as it determines whichtransform of the model parameter should be modeled as linear in the observed each of the three examples discussed, we used what is called thenatural link, whichis motivated by considering a change-of-variable for the parameter, 7 ( ), so that thePDF/PMFf(y| ) in terms of the new parameter has the formf(y| ) =e y A( )h(y)for some functionsAandh.

9 For example, the Bernoulli PMF isf(y) =py(1 p)1 y= (1 p)(p1 p)y=e(logp1 p)y+log(1 p),so we may set = logp1 p,A( ) = log(1 p) = log(1 +e ), andh(y) = 1. This is calledtheexponential familyform of the PDF/PMF, and is called thenatural each example, the natural link simply setsg( ) = ( ) (or equivalently,g( ) =c ( ) fora constantc).Use of the natural link leads to some nice mathematical properties for likelihood-basedinference for instance, since is modeled as linear in , the second-order partial derivativesoflogf(Y| ) = Y A( ) + logh(Y)27-3with respect to do not depend onY, so the Fisher information is always given by 2l( )without needing to take an expectation.

10 (We sometimes say in this case that the observedand expected Fisher information matrices are the same.) On the other hand, from themodeling perspective, there is usually no intrinsic reason to believe that the natural linkg( ) = ( ) is the correct transformation of that is well-modeled as a linear combinationof the covariates, and other link functions are also commonly used, especially if they lead toa better fit for the