# Parameter estimation in a drifting normal distribution

Suppose we observe $$n$$ pairs of points $$(a_1,b_1),~(a_2,b_3)~..~(a_n,b_n)$$. The underlying data generating process is known to be as follows:

• $$u_i \sim N(0,1)$$ and $$a_i \sim N(u_i,1)$$
• $$b_i \sim N(u_i + k*a_i, 1)$$

independently for each $$i$$ and $$k$$ is a constant. What is a good unbiased estimator for $$k$$?

• We have $u = n_0$, $a = n_0 + n_1$, $b = n_0 + k*(n_0 + n_1)+ n_2 = (1+k)*n_0 + k*n_1 + n_2$ where $n_i$'s are all standard normals. Variance of $b$ then is $(1+k)^2+k^2+1$. Equating this with empirical estimate of variance and solving gives you a pretty good estimate for $k$. But it is not unbiased.
– user284807
Jul 25, 2022 at 23:13
• @bleh So simple, but so efficient... I developped formulas based on quadratic optimisations of bayesian estimators...
– FP0
Jul 25, 2022 at 23:19
• Haha. It would be interesting to see how to estimate the drift in general. Say if we have $b \sim N(u+f(a), 1)$ for smooth function $f$ and we need to estimate it. Maybe your methods will come handy then?
– user284807
Jul 25, 2022 at 23:25
• I tried to write a very comprehensive answer. @bleh, user363943, may I kindly ask you to help me to check that what I wrote is correct and contains no typo?
– FP0
Jul 26, 2022 at 1:31

## 1 Answer

I am going to start to discuss @bleh's idea in order to explain why my solution seems more robust.

1. @bleh's idea:

Assume:

• $$N_{i,j} \sim N(\mu_j,\sigma^2_j)$$ are independent random variables $$\forall i,j$$
• $$U_i=N_{i,1}$$
• $$A_i\sim N(U_i,\sigma_2^2)\Rightarrow A_i=N_{i,1}+N_{i,2}$$
• $$B_i \sim N(U_i+kA_i,\sigma_3^2)\Rightarrow B_i=(1+k)N_{i,1}+kN_{i,2}+N_{i,3}$$

As a consequence, $$B_i\sim N((1+k)\mu_1+k\mu_2+\mu_3,k^2(\sigma^2_1+\sigma^2_2)+2k\sigma^2_1+\sigma^3_1)$$

It is tempting to compute an estimator $$\widehat{\sigma^2_B}$$ of the variance of $$B$$ based on the squared differences from the means $$\sum_{i=1}^n\left(b_i-\mathbb{E}(B_i)\right)^2$$, to equate it to $$\mathbb{V}ar(B)$$, and to solve for $$k$$.

This is a problem for 2 main reasons:

1. To solve this polynomial, the value of the discriminant $$\Delta$$ will depend on the parameters $$\sigma^2_j$$ as well as on the estimator $$\widehat{\sigma^2_B}$$. However, in some unlikely cases (ex: very small $$\widehat{\sigma^2_B}$$ being computed despite $$\sigma^2_3$$ being large), the calculation will fail (due to $$\Delta<0$$) because the observed values $$b_i$$ were not dispersed enough compared to the expected variance $$\sigma^2_3$$.
2. Even if this process worked, determining an estimator $$\widehat{\sigma^2_B}$$ that would allow us to find an unbiased estimator for $$k$$ after all this steps, including the calculations of squared deviations to a mean, and a square root, would be very difficult.

In addition, we would end up loosing some information because we would use our observations $$b_i$$, but not $$a_i$$.

Hence, I propose my solution:

1. Use bayesian statistics and the normal conjugate prior.

You can see here that:

If $$\mu \sim N(\mu_0,\sigma^2_0)$$ and $$X_i\sim N(\mu,\sigma^2)$$ (the expected value of $$X_i$$ depends on a random value $$\mu$$), then:

• $$\mu\vert x_i, i=1,...,n \sim N(\mu_0\prime, \sigma^{2}_0\prime)$$, with $$\mu_0\prime = \frac{1}{\frac{1}{\sigma^2_0}+\frac{n}{\sigma^2}}\left(\frac{\mu_0}{\sigma^2_0}+\frac{\sum_{i=1}^{n}x_i}{\sigma^2}\right)$$ and $$\sigma_0^2\prime=\frac{1}{\frac{1}{\sigma_0^2}+\frac{n}{\sigma^2}}$$ (more observations of $$x_i$$ allow to estimate a more accurate distribution for $$\mu$$)
• $$X_{i+1} \sim N(\mu_0\prime, \sigma^{2}_0\prime + \sigma^2)$$ (updated knowledge about $$\mu$$ allow to estimate a more accurate distribution for future observations of $$X$$)

Apply this with your example:

• $$U_i\sim N(\mu_U,\sigma_U^2)$$
• $$A_i \sim N(U_i, \sigma_A^2)$$
• $$B_i \sim N(U_i+kA_i, \sigma_B^2)$$

So this is a "tower" of normal variables using the previous variables as conjugate prior.

The observations we have are the couples $$\{a_i,b_i\}$$.

Let's start by the end:

$$B_i\vert A_i \sim N(U_i+kA_i, \sigma_B^2)\vert A_i$$

In this "equation", we know $$A_i$$, but not $$U_i$$. But based on our bayesian knowledge, we can determine the distribution of $$U_i \vert A_i$$:

Based on the equations above:

$$U_i \lvert A_i \sim N\left(\tilde{\mu_i}=\frac{1}{\frac{1}{\sigma_U^2}+\frac{1}{\sigma_A^2}}\left(\frac{\mu_U}{\sigma_U^2}+\frac{A_i}{\sigma_A^2}\right),\tilde{\sigma^2_{U,i}}=\frac{1}{\frac{1}{\sigma_U^2}+\frac{1}{\sigma^2_A}}\right)$$

Hence, you will have:

$$B_i\vert A_i \sim N(kA_i+\tilde{\mu_i},\sigma^2_B+\tilde{\sigma^2_{U,i}})$$

Hence, $$\mathbb{E}(B_i\vert A_i)=kA_i+\tilde{\mu_i}$$, and we can propose a range of estimators $$\hat{k}_i=\frac{b_i-\tilde{\mu_i}}{a_i}$$. They are unbiased, and each has a variance which can be computed with the parameters above.

Since you have $$i$$ estimators $$\hat{k}_i$$ with variance $$\sigma^2_{k, i}=\frac{\sigma^2_B+\tilde{\sigma^2_{U,i}}}{a_i^2}$$, you can determine a linear combination of these estimators minimising the variance of the final estimator:

$$\min_w w' \Sigma w$$, st $$w'\mathbb{1}=1$$ where $$\Sigma$$ is the covariance matrix of the estimators $$\hat{k}_i$$.

Finally, the unbiased estimator with the smallest variance will be:

$$w^*=\frac{\Sigma^{-1}\mathbb{1}}{\mathbb{1}'\Sigma^{-1}\mathbb{1}}$$