This is my first Stack Exchange post - I hope I've followed the guidelines correctly.
There have been a few previous questions asked about negative variance, from which I've learned some things, but none of them have quite addressed my issue.
In the image below is some data I generated via simulation (stripped of all physics to gain acceptance here :P )
I need to use this simulation to analyse a separate set of data, which are measurements. I don't need the fit of the simulation data to be physical, but just analytical so I can interpolate. As such, I used Adj. $R^{2}$ tests to settle on a 4th-order polynomial. It's not perfect, but it's good enough.
As a result of the fitting parameters all being highly correlated, I've attempted to use the full matrix approach to propagate uncertainty. The variance-covariance matrix was generated by the fitting software (Origin), and the nature of a 4th-order polynomial makes the Jacobian quite simple:
$$ y = B_{0} +B_{1}x + B_{2}x^{2} + B_{3}x^{3} + B_{4}x^{4} \mathrm{,and} $$
$$ \Sigma^{y} = \mathbf{J\Sigma^{B}J}^{\mathrm{T}} \mathrm{, where}$$
$$ \mathbf{J} = \left[ \begin{array}{cc} 1 & x & x^2 & x^3 & x^4 \end{array} \right] \mathrm{, and}$$
$$ \mathbf{\Sigma^{B}} = \left[ \begin{array}{cc} 0.00413 & -0.000185 & 2.95\times 10^{-6} & -2.01\times 10^{-8} & 4.93\times 10^{-11} \\ -0.000185 & 8.31\times 10^{-6} & -1.34\times 10^{-7} & 9.13\times 10^{-10} & -2.25\times 10^{-12} \\ 2.95\times 10^{-6} & -1.34\times 10^{-7} & 2.16\times 10^{-9} & -1.48\times 10^{-11} & 3.68\times 10^{-14} \\ -2.01\times 10^{-8} & 9.13\times 10^{-10} & -1.48\times 10^{-11} & 1.02\times 10^{-13} & -2.55\times 10^{-16} \\ 4.93\times 10^{-11} & -2.25\times 10^{-12} & 3.68\times 10^{-14} & -2.55\times 10^{-16} & 6.38\times 10^{-19} \end{array} \right] \mathrm{.} $$
I know that, under normal circumstances, many of the terms in the above would be so low as to be negligible. However, the terms in the correlation matrix are all >0.9, and the size of the terms in the Jacobian cancel out much of the "negligibility". I understand that it's fine for the matrix elements themselves to be negative, but that $\Sigma^{y}$ should not be.
Working this through for all $x$ gives negative values for the variance at each point. For example, when $x = 50$, $\Sigma^{y} = -0.00013$. This is problematic in itself, but particularly as I was hoping to use a real-valued standard deviation as my error bars!
Any help you can give me in untangling this is greatly appreciated. Normally, I deal with uncorrelated variables, so this is a bit out of my comfort zone. I am specifically looking to know:
- If you believe my approach is correct,
- If I am missing something trivial, e.g. whether variance can actually be negative in this situation, and
- Any advice on how to proceed.
Thank you in advance!
So, as user619894 posted in the comments, the problem was numerical. I had copied the covariance matrix into excel for processing, and in doing so had truncated the values. Changing the format of the cells and re-pasting such that a larger number of sig. figs. was used, fixed the issue and turned the variance positive.