Point estimation

Estimation methods

Below are some commonly used methods for estimating unknown parameters. The methods vary in their domain of applicability and on the underlying statistical paradigm (frequentist or Bayesian).

Frequentist methods

Least square estimation (LSE)

The method of least squares is a parameter estimation procedure popular for regression, where observations ${\displaystyle Y_1,\dots ,Y_n}$ need to be expressed as functions of the covariates ${\displaystyle X_1,\dots ,X_n\in {\mathbb{R}}^p}$ via

$${\displaystyle Y_i=f_{\beta }(X_i)+{\epsilon }_i}$$

for some function ${\displaystyle f_{\beta }}$ depending on a parameter vector ${\displaystyle \beta \in {\mathbb{R}}^{p+1}}$ and a random variable ${\displaystyle {\epsilon }_i}$ accounting for random noise in the data. The ordinary least squares estimator (OLS) ${\displaystyle {\hat{\beta }}_{\mathrm{O}\mathrm{L}\mathrm{S}}}$ assumes that ${\displaystyle f_{\beta }}$ is a linear function of ${\displaystyle \beta }$ and estimates the parameter by minimizing the sum of squared differences, called residuals:

$${\displaystyle {\hat{\beta }}_{\mathrm{O}\mathrm{L}\mathrm{S}}\in \mathrm{a}\mathrm{r}\mathrm{g}\underset{\beta \in {\mathbb{R}}^p}{min}\sum _{i=1}^n{\left(Y_i-{\beta }_0-X_{i,1}{\beta }_1-\cdots -X_{i,p}{\beta }_p\right)}^2.}$$

In fact, if the noise variables ${\displaystyle {\epsilon }_i}$ are independent and identically normally distributed with mean zero, then the OLS coincides with the maximum likelihood estimator.

If the rank of the matrix of covariates ${\displaystyle X=(X_1,\dots ,X_n{)}^{\mathrm{\top }}\in {\mathbb{R}}^{n\times p}}$ is equal to ${\displaystyle p\le n}$, then the OLS is uniquely given by

$${\displaystyle {\hat{\beta }}_{\mathrm{O}\mathrm{L}\mathrm{S}}=(X^{\mathrm{\top }}X{)}^{-1}{X}^{\mathrm{\top }}Y,}$$

and it is the best linear unbiased estimator (BLUE) by the Gauß–Markov theorem. If, as above, the noise variables are independent and identically normally distributed, then the OLS is even UMVU for the entire class of unbiased estimators.

There are multiple variants of least squares estimation:

• weighted least squares (WLS) estimation: The residuals are weighted with different weighting factors ${\displaystyle w_i}$.
• generalized least squares (GLS) estimation: The noise variables are allowed to be correlated.
• nonlinear least squares estimation: The function ${\displaystyle f_{\beta }}$ can depend non-linearly on the parameter ${\displaystyle \beta }$.
• regularized least squares estimation: The parameter is chosen such that the sum of squared residuals plus a regularization term is minimized. Typical choices for the regularization are the euclidian norm ${\displaystyle \Vert \beta {\Vert }_2^2}$, leading to Ridge regression, or the 1-norm ${\displaystyle \Vert \beta {\Vert }_1}$, leading to Lasso regression.

Maximum likelihood estimation (MLE)

The method of maximum likelihood, due to R.A. Fisher, is arguably the most important general method for estimating the parameters ${\displaystyle {\theta }_1,\dots ,{\theta }_p}$ of a parametric model. According to this method, the best estimate for the unknown model parameter is the one that maximizes the probability (or probability density) of observing the data that has been observed. This probability as a function of the model parameters is called the likelihood, giving the method its name.

In mathematical terms, the methods works as follows: Let ${\displaystyle X=(X_1,\dots ,X_n)}$ denote a random data sample with joint probability density function or probability mass function ${\displaystyle f_{\theta }}$ depending on the vector ${\displaystyle \theta =({\theta }_1,\dots ,{\theta }_p)}$ of model parameters. The function ${\displaystyle \theta \mapsto f_{\theta }(X)}$ is called the likelihood function, often denoted by ${\displaystyle L}$. Assuming that true model parameter lies in a set ${\displaystyle \mathrm{\Theta }\subset {\mathbb{R}}^p}$, a maximum likelihood estimator ${\displaystyle {\hat{\theta }}_{\mathrm{M}\mathrm{L}\mathrm{E}}}$ fulfills

$${\displaystyle L({\hat{\theta }}_{\mathrm{M}\mathrm{L}\mathrm{E}})=\underset{\theta \in \mathrm{\Theta }}{max}L(\theta )=\underset{\theta \in \mathrm{\Theta }}{max}{f}_{\theta }(X).}$$

In practice, the likelihood function is often differentiable, in which case ${\displaystyle {\hat{\theta }}_{\mathrm{M}\mathrm{L}\mathrm{E}}}$ is a solution of the equations

$${\displaystyle \frac{\mathrm{\partial }L(\theta )}{\mathrm{\partial }{\theta }_i}=0,i=1,\dots ,k.}$$

Under certain assumptions on the likelihood, the MLE is strongly consistent, asymptotically efficient and asymptotically normal.

Method of moments (MoM)

The method of moments is one of the oldest methods of estimation, which was first introduced by Karl Pearson in 1895. This method uses the functional relationship between quantities of interest of a distribution (e.g. the variance or model parameters) and the moments of that distribution, leading to a set of equations that are solved for the aforementioned quantities. Motivated by the law of large numbers, the moments of the distribution are estimated by the sample moments, that is, the arithmetic averages of powers of the data values. The method stands out for its simplicity and wide domain of applicatibility, but often do not yield the best estimators and sometimes even no estimator at all. There also exists generalizations of the method of moments.

Mathematically speaking, for a distribution ${\displaystyle {\mathbb{P}}_{\theta }}$ depending on some vector of parameters ${\displaystyle \theta =({\theta }_1,\dots ,{\theta }_k)}$, the equations of the first ${\displaystyle k}$ moments of ${\displaystyle X\sim {\mathbb{P}}_{\theta }}$ are established:

$${\displaystyle \begin{array}{cc}& m_1=\mathbb{E}[X]=f_1({\theta }_1,\dots ,{\theta }_k),\\ & \dots \\ & m_k=\mathbb{E}[X^k]=f_k({\theta }_1,\dots ,{\theta }_k).\end{array}}$$

If these equations have a solution, then there exists functions ${\displaystyle g_1,\dots ,g_k}$ such that the parameters can be expressed by the moments via

$${\displaystyle \begin{array}{cc}& {\theta }_1=g_1(m_1,\dots ,m_k),\\ & \dots \\ & {\theta }_k=g_k(m_1,\dots ,m_k).\end{array}}$$

Given a random sample ${\displaystyle X_1,\dots ,X_n}$ from ${\displaystyle {\mathbb{P}}_{\theta }}$, the moments are then replaced by the sample moments ${\displaystyle {\hat{m}}_j=\frac{1}{n}\sum _{i=1}^n{X}_i^j}$, ${\displaystyle j=1,\dots ,k}$, yielding the moment estimators:

$${\displaystyle \begin{array}{cc}& {\hat{\theta }}_1=g_1({\hat{m}}_1,\dots ,{\hat{m}}_k),\\ & \dots \\ & {\hat{\theta }}_k=g_k({\hat{m}}_1,\dots ,{\hat{m}}_k).\end{array}}$$

Under relatively low assumptions on the solution functions ${\displaystyle g_1,\dots ,g_k}$ and the moments of ${\displaystyle {\mathbb{P}}_{\theta }}$, the moment estimators are strongly consistent and asymptotically normal.

Lehmann–Scheffé theorem

The Lehmann–Scheffé theorem, although not an estimation method in itself, motivates two methods for constructing UMVU estimators in models, where a complete, sufficient statistic ${\displaystyle T}$ is available. This applies, for example, to exponential families, where the natural sufficient statistic is always complete.

Let ${\displaystyle X=(X_1,\dots ,X_n)}$ be a vector of random samples from some distribution ${\displaystyle {\mathbb{P}}_{\theta }}$ depending on a parameter ${\displaystyle \theta }$.

Method 1: For the estimation of ${\displaystyle \theta }$, there exists at most one function ${\displaystyle g}$ that solves the equation $${\displaystyle {\mathbb{E}}_{\theta }[g(T(X))]=\theta ,}$$ in which case the estimator ${\displaystyle \hat{\theta }=g(T)}$ is an UMVU estimator for ${\displaystyle \theta }$.

Method 2: Given any unbiased estimator ${\displaystyle \tilde{\theta }}$ of ${\displaystyle \theta }$, the estimator $${\displaystyle \hat{\theta }(X)={\mathbb{E}}_{\theta }[\tilde{\theta }(X)\mid T(X)],}$$ that is, the conditional expectation of ${\displaystyle \tilde{\theta }}$ given ${\displaystyle T}$, is the unique UMVU estimator. In some cases, computing the conditional expectation of a suitably chosen unbiased estimator given ${\displaystyle T}$ might be easier than solving the equation from method 1.

Score matching estimation

Score matching is a rather new method for parameter estimation that is less intuitive than MLE, but is more advantageous for complicated distributions with a large number of parameters.

The data is assumed to be drawn from a distribution with probability density ${\displaystyle f(\mathbf{x},\mathit{\theta })}$, ${\displaystyle \mathbf{x}\in {\mathbb{R}}^d}$, parameterized by ${\displaystyle \mathit{\theta }}$. The key quantity, that is needed for the score matching estimation, is the score ${\displaystyle {\mathrm{\nabla }}_{\mathbf{x}}\log f(\mathbf{x},\mathit{\theta })}$. The true parameter (vector) ${\displaystyle \mathit{\theta }}$ is characterized by$${\displaystyle \begin{array}{cc}\mathit{\theta }& =\mathrm{a}\mathrm{r}\mathrm{g}\underset{{\mathit{\theta }}^{\prime }}{min}\frac{1}{2}{\int }_{{\mathbb{R}}^d}\Vert {\mathrm{\nabla }}_{\mathbf{y}}\log f(\mathbf{y},{\mathit{\theta }}^{\prime })-{\mathrm{\nabla }}_{\mathbf{y}}\log f(\mathbf{y},\mathit{\theta }){\Vert }^{2}f(\mathbf{y},\mathit{\theta })\mathrm{d}\mathbf{y}\\ & =\mathrm{a}\mathrm{r}\mathrm{g}\underset{{\mathit{\theta }}^{\prime }}{min}\sum _{i=1}^d{\int }_{{\mathbb{R}}^d}\left({\mathrm{\partial }}_{y_i}^2\log f(\mathbf{y},{\mathit{\theta }}^{\prime })+\frac{1}{2}({\mathrm{\partial }}_{y_i}\log f(\mathbf{y},{\mathit{\theta }}^{\prime }{)}^2\right)f(\mathbf{y},\mathit{\theta })\mathrm{d}\mathbf{y}.\end{array}}$$In line with the second equality, for observations ${\displaystyle {\mathbf{x}}_1,\dots ,{\mathbf{x}}_n}$, the score matching estimate ${\displaystyle \hat{\mathit{\theta }}}$ is defined by$${\displaystyle \hat{\mathit{\theta }}=\mathrm{a}\mathrm{r}\mathrm{g}\underset{{\mathit{\theta }}^{\prime }}{min}\sum _{j=1}^n\sum _{i=1}^d\left({\mathrm{\partial }}_{x_i}^2\log f({\mathbf{x}}_j,{\mathit{\theta }}^{\prime })+\frac{1}{2}({\mathrm{\partial }}_{x_i}\log f({\mathbf{x}}_j,{\mathit{\theta }}^{\prime }{)}^2\right).}$$The main advantage of score matching estimation compared to the maximum likelihood estimation is the following. Typically, the probability density is of the form$${\displaystyle f(\mathbf{x},\mathit{\theta })=\frac{1}{Z(\mathit{\theta })}q(\mathbf{x},\mathit{\theta }),}$$with an unnormalized density ${\displaystyle q(\mathbf{x},\mathit{\theta })}$ and a normalization constant ${\displaystyle Z(\mathit{\theta })}$. The computation of latter is often the bottleneck of maximum likelihood estimation, especially if the dependence on ${\displaystyle \mathit{\theta }}$ is non-analytic and the parameter vector is high-dimensional. By taking the derivatives of the logarithm of ${\displaystyle f}$, ${\displaystyle Z(\mathit{\theta })}$ does not occur in the score matching estimation procedure, thus eliminating this issue.

Additionally, under mild restrictions on ${\displaystyle f(\mathbf{x},\mathit{\theta })}$, the score matching estimate can be shown to be weakly consistent.

Bayesian methods

Bayesian estimation methods take into account the statistician's prior belief on the distribution of the parameters, which is modeled as a distribution with density ${\displaystyle \pi }$ on the parameter space ${\displaystyle \mathrm{\Theta }}$, called the prior distribution. Given observations ${\displaystyle X_1,\dots ,X_n}$, which, for fixed ${\displaystyle \theta \in \mathrm{\Theta }}$, are assumed to be distributed according to a distribution with density ${\displaystyle p_{\theta }}$, the prior distribution is updated according to Bayes' rule:

$${\displaystyle p(\theta \mid X_1,\dots ,X_n)=\frac{p_{\theta }(X_1,\dots ,X_n)\pi (\theta )}{{\int }_{\mathrm{\Theta }}{p}_{\theta }(X_1,\dots ,X_n)\pi (\theta )\mathrm{d}\theta }.}$$The updated distribution of the parameters given the data is called the posterior distribution. Many Bayesian point estimators are the posterior distribution's statistics of central tendency, e.g., its mean, median, or mode.

Bayes estimator

Given a vector of random observations ${\displaystyle X=(X_1,\dots ,X_n)}$ and a loss function ${\displaystyle L}$, a Bayes estimator is any estimator ${\displaystyle \hat{\theta }(X)}$ that minimizes the Bayes risk: $${\displaystyle R_{\pi }(\hat{\theta })={\int }_{\mathrm{\Theta }}\mathbb{E}[L(\theta ,\hat{\theta }(X))]\pi (\theta )\mathrm{d}\theta =\mathbb{E}\left[{\int }_{\mathrm{\Theta }}L(\theta ,\hat{\theta }(X))p(\theta \mid X)\mathrm{d}\theta \right].}$$

Obviously, the risk is minimized when the integral under the expectation is minimized as a function of the observations. Thus, a Bayes estimator is equivalently given as a minimizer of the posterior loss:

$${\displaystyle \hat{\theta }(X)\in \mathrm{a}\mathrm{r}\mathrm{g}\underset{y\in {\mathbb{R}}^n}{min}{\int }_{\mathrm{\Theta }}L(\theta ,y)p(\theta \mid X)\mathrm{d}\theta .}$$

The most common choices for ${\displaystyle L}$ are the following:

• ${\displaystyle L(\theta ,y)=(\theta -y{)}^2}$: In this case, the posterior loss is minimized by the posterior mean ${\displaystyle \hat{\theta }(X)={\int }_{\mathrm{\Theta }}\theta p(\theta \mid X)\mathrm{d}\theta }$, as observed by Gauss.

• ${\displaystyle L(\theta ,y)=|\theta -y|}$: Here, the posterior loss is minimized by the median of the posterior distribution, as observed by Laplace.

By Wald's theorem, unique Bayesian estimators are admissible. Moreover, under certain assumptions on the prior and ${\displaystyle p_{\theta }}$, the posterior mean estimator can be shown to be consistent, asymptotically normal and asymptotically efficient.

Maximum a posteriori (MAP)

The maximum a posteriori estimator ${\displaystyle {\hat{\theta }}_{\mathrm{M}\mathrm{A}\mathrm{P}}}$ is a maximizer of the posterior distribution:

$${\displaystyle {\hat{\theta }}_{\mathrm{M}\mathrm{A}\mathrm{P}}(X)\in \mathrm{a}\mathrm{r}\mathrm{g}\underset{\theta \in \mathrm{\Theta }}{max}p(\theta \mid X).}$$

In a sense, the MAP estimator is a Bayesian analog to the maximum likelihood estimator. In fact, for a uniform prior distribution, the MAP estimator coincides with the MLE.

The MAP estimator has good asymptotic properties, even for many difficult problems, on which the maximum-likelihood estimator has difficulties. For regular problems, where the maximum-likelihood estimator is consistent, the maximum-likelihood estimator ultimately agrees with the MAP estimator.

Others

The Minimum Message Length (MML) point estimator is based in Bayesian information theory and is not so directly related to the posterior distribution.

Special cases of Bayesian filters are important:

• Kalman filter
• Wiener filter

Several methods of computational statistics have close connections with Bayesian analysis:

• particle filter
• Markov chain Monte Carlo (MCMC)

Point estimate v.s. confidence interval estimate

There are two major types of estimates: point estimate and confidence interval estimate. In the point estimate we try to choose a unique point in the parameter space which can reasonably be considered as the true value of the parameter. On the other hand, instead of unique estimate of the parameter, we are interested in constructing a family of sets that contain the true (unknown) parameter value with a specified probability. In many problems of statistical inference we are not interested only in estimating the parameter or testing some hypothesis concerning the parameter, we also want to get a lower or an upper bound or both, for the real-valued parameter. To do this, we need to construct a confidence interval.

Confidence interval describes how reliable an estimate is. We can calculate the upper and lower confidence limits of the intervals from the observed data. Suppose a dataset x₁, . . . , x_n is given, modeled as realization of random variables X₁, . . . , X_n. Let θ be the parameter of interest, and γ a number between 0 and 1. If there exist sample statistics L_n = g(X₁, . . . , X_n) and U_n = h(X₁, . . . , X_n) such that P(L_n < θ < U_n) = γ for every value of θ, then (l_n, u_n), where l_n = g(x₁, . . . , x_n) and u_n = h(x₁, . . . , x_n), is called a 100γ% confidence interval for θ. The number γ is called the confidence level. In general, with a normally-distributed sample mean, Ẋ, and with a known value for the standard deviation, σ, a 100(1-α)% confidence interval for the true μ is formed by taking Ẋ ± e, with e = z_1-α/2(σ/n^1/2), where z_1-α/2 is the 100(1-α/2)% cumulative value of the standard normal curve, and n is the number of data values in that column. For example, z_1-α/2 equals 1.96 for 95% confidence.

Here two limits are computed from the set of observations, say l_n and u_n and it is claimed with a certain degree of confidence (measured in probabilistic terms) that the true value of γ lies between l_n and u_n. Thus we get an interval (l_n and u_n) which we expect would include the true value of γ(θ). So this type of estimation is called confidence interval estimation. This estimation provides a range of values which the parameter is expected to lie. It generally gives more information than point estimates and are preferred when making inferences. In some way, we can say that point estimation is the opposite of interval estimation.