I want to calculate numerically the maximum likelihood estimators of $(\beta,\sigma)$ for the linear regression model:

$$y_j = x_j^{\top}\beta + \epsilon_j, $$

where $j=1,\dots,n$, $\beta$ is $p$-dimensional, and $\epsilon_j$ are i.i.d. according to a Laplace distribution with location zero and scale $\sigma$.

Given the differentiability issues with this distribution, I cannot obtain the MLE of $\beta$ by solving the score functions. I need to do this for a number of data sets with large $p$. Is there an efficient method for doing so?


The MLE minimizes

$\min \| X \beta - y\|_{1}$

Unfortunately, there isn't any simple closed form solution to this optimization problem. However, this is a convex optimization problem for which there are many available approaches.

For reasonably small instances (e.g. $n$ on the order of $100,000$ and $p$ is less than say $1,000$), the simplest approach is to use a linear programming library routine to solve the linear programming problem:

$\min \sum_{i=1}^{n} t_{i}$

subject to

$t \geq X\beta-y$

$t \geq y-X\beta$.

For larger instances (e.g. $n$ is in the billions or even larger), you might consider a stochastic subgradient descent method. This is relatively easy to implement in a "big data" environment with hadoop.

  • $\begingroup$ Many thanks for your answer. Your suggestions look very interesting. $\endgroup$ – TheRock Apr 8 '15 at 14:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.