# Why isn't the pseudo-inverse the best choice in my linear estimation problem?

Context:

I have a problem of the following form. Let $$\boldsymbol{\theta}\in\mathbb{R}^n$$ be a fixed vector I want to estimate. Let $$\mathbf{M}\in\mathbb{R}^{m\times n}$$ be a matrix with $$m>n$$ and $$\text{rank}(\mathbf{M})=n$$. Moreover, let $$\mathbf{w}\in\mathbb{R}^m$$ be a vector with independent identically distributed components, each distributed as $$\sim\mathcal{N}(0,\sigma^2)$$. Assume that you have access to the "noisy vector": $$\mathbf{x} = \mathbf{M}\boldsymbol{\theta} + \mathbf{w}$$ The question is to obtain the "best" estimate $$\hat{\boldsymbol{\theta}}$$ of $$\boldsymbol{\theta}$$, under mean squared error criterion. My first instinct is to set my estimate using a pseudo-inverse according to $$\hat{\boldsymbol{\theta}}_1 = (\mathbf{M}^\top\mathbf{M} )^{-1}\mathbf{M}^\top\mathbf{x}$$

However, in a particular problem I have, I found by "accident" another matrix $$\mathbf{M}_{\text{special}}$$ which seems to perform better than $$(\mathbf{M}^\top\mathbf{M} )^{-1}\mathbf{M}^\top$$ somehow. Call $$\hat{\boldsymbol{\theta}}_2 = \mathbf{M}_{\text{special}}\mathbf{x}$$ for this special matrix I found. The concrete data I'm using can be found by running the matlab script at the end of this question.

To test which estimate performs better, I'm computing the "overall" errors $$\sqrt{\frac{1}{n}({\boldsymbol{\theta}}-\hat{\boldsymbol{\theta}}_1)^\top({\boldsymbol{\theta}}-\hat{\boldsymbol{\theta}}_1)}$$ and $$\sqrt{\frac{1}{n}({\boldsymbol{\theta}}-\hat{\boldsymbol{\theta}}_2)^\top({\boldsymbol{\theta}}-\hat{\boldsymbol{\theta}}_2)}$$ and comparing them. My observations are the following:

• If $$\sigma^2=0$$ (no noise), then I can recover $$\hat{\boldsymbol{\theta}}_1={\boldsymbol{\theta}}$$, which is expected since $$\mathbf{M}^\top\mathbf{M}$$ is full rank. On the other hand $${\boldsymbol{\theta}}-\hat{\boldsymbol{\theta}}_2$$ is a non-zero vector (even without noise). Hence, the pseudo-inverse estimate performs better.
• As I increase the noise variance, both estimates decrease performance (which is expected). However, the pseudo-inverse one decreases faster. In fact, when $$\sigma=0.5$$, both estimates perform roughly the same.
• For bigger values of $$\sigma^2>1$$, I find that $$\hat{\boldsymbol{\theta}}_2$$ ALWAYS performs better than $$\hat{\boldsymbol{\theta}}_1$$.

Question:

How is it possible for a different linear estimate $$\hat{\boldsymbol{\theta}}_2$$ outperform what is supposed to be the best linear estimator $$\hat{\boldsymbol{\theta}}_1$$? I assume I misunderstood something about the "optimality nature" of the pseudo-inverse estimate I'm using. However, I cannot find the error in my reasoning. What makes $$\hat{\boldsymbol{\theta}}_2$$ better?

Data:

The data and script I'm using can be found here: https://github.com/RodrigoAldana/se_question_data

I know that I'm generating new $$\boldsymbol{\theta}$$ in each run. Please note that I also tried to fix $$\boldsymbol{\theta}$$ for all runs obtaining similar results.

• Since $\theta$ is a vector, how are you calculating the standard deviation of $\theta - \hat{\theta}$? Commented Jan 16, 2023 at 22:57
• And why do you need a pseudo-inverse? Your design matrix $M$ is of full rank, by assumption, meaning there is a unique inverse of $M^\prime M$ and it is equal to any pseudo-inverse.
– whuber
Commented Jan 16, 2023 at 22:59
• $(M'M)^{-1}M'x$ is the best linear unbiased estimator, and has a covariance matrix $\sigma^2(M'M)^{-1}$, which differs from all other unbiased linear estimators by a positive semidefinite matrix. Your estimator is evidently biased, among other things. Without knowing how your "special matrix" was generated, it's extremely difficult (for me, at any rate) to say more. It is well-known that biased linear estimators, e.g., Ridge regression estimators, can outperform OLS. Commented Jan 16, 2023 at 23:07
• The "U" in BLUE stands for unbiased. It is well-known that biased linear estimators, e.g., Ridge regression estimators, can outperform OLS. There is no contradiction. Commented Jan 16, 2023 at 23:11
• No worries! We all have simple questions sometimes! Commented Jan 16, 2023 at 23:14

$$(M'M)^{-1}M'x$$ is the best linear unbiased estimator and has a covariance matrix $$\sigma^2(M′M)^{−1}$$, which differs from all other unbiased linear estimators by a positive semidefinite matrix. Your estimator is biased, possibly among other things, so the Gauss-Markov theorem does not apply to it.