Consider a time evolving phenomenon represented by a variable $y(t)$, whose dynamics is dependent on a parameter, say the temperature $\theta$.

We have two series of measurements at different constant temperatures $\theta_1$ and $\theta_2$. The plot of look like this:

Sketch of plot.

Measurements at a given temperature cannot be considered independent because of say the imprecision of the control over experimental conditions.

Measurements for different temperatures can reasonably be considered independent.

We would like to be able to predict the time evolution of $y$ at different temperatures.

We model the dynamics of $y$ as $y(t) = a(\theta) t$, to be used later on as $\dot{y} = a(\theta)$.

Assuming an additive error term, the estimation of $a$ is a standard regression problem that can be solved by least square optimisation. The main issue is the limited amount of data, especially when considering the intra-series dependency. Thankfully, some physical consideration tell us that $a$ varies exponentially with the temperature: $a(\theta) = \alpha \: e^{\beta / \theta}$.

Now we are interested into modeling the uncertainty on $\hat{\alpha}$ and $\hat{\beta}$ the estimators of $\alpha$ and $\beta$. What procedure would you recommend?

One idea would be to use bootstrap after non-linear least square fitting. However, we are mainly interested into the extrapolation area, namely predicting the behaviour of $y$ for $\theta \notin [\theta_1, \theta_2]$ (actually $\theta \leq \theta_1 \leq \theta_2$, and we could probably provide a lower bound for $a$).

  • $\begingroup$ Unfortunately your link is not working for me. $\endgroup$ – Glen_b -Reinstate Monica Jun 10 '16 at 8:04
  • $\begingroup$ @Glen_b: maybe the picture insertion feature of stackexchange require a little delay? $\endgroup$ – Alfred M. Jun 10 '16 at 8:19
  • $\begingroup$ I am getting the message that the image contains errors and cannot be displayed. $\endgroup$ – Glen_b -Reinstate Monica Jun 10 '16 at 8:49

Here is where I would start.

Bayesian inference for non-linear regression model

Let $y_{i,t}$ be the observation at time $t$ with temperature $\theta_i$ where $i=1,2$ (since there are two temperatures) and $t=1,\ldots,T$ (for simplicity I'm assuming observations are taken at the same times, but this could easily be modified for situations where this is not the case.) Assume an additive normal error, i.e.

$$ y_{i,t} \stackrel{ind}{\sim} N(\alpha e^{\beta/\theta_i}, \sigma^2). $$

A Bayesian approach to estimation requires a prior over the parameters, i.e. $p(\alpha,\beta,\sigma^2)$, and then the posterior is

$$ p(\alpha,\beta,\sigma^2|y) \propto p(\alpha,\beta,\sigma^2) \prod_{i=1}^2 \prod_{t=1}^{T} N\left(y_{i,t}; \alpha e^{\beta/\theta_i} t, \sigma^2\right). $$ which will likely need to be estimated computationally, e.g. Markov chain Monte Carlo.

This posterior can provide point estimates, e.g. $E[\alpha|y]$, and uncertainties, e.g. $V[\alpha|y]$. Forecasts can be obtained for a new temperature $\tilde{\theta}$ and time $\tilde{t}$. Assuming the response is independent of the previous data given the parameters, the forecast distribution is

$$ p(\tilde{y}|y) = \int \int \int N\left(\tilde{y};\alpha e^{\beta/\tilde{\theta}} \tilde{t}, \sigma^2\right) p(\alpha,\beta,\sigma^2|y) d\alpha d\beta d\sigma^2 $$

which will also likely need to be estimated computationally.

Although this forecast distribution will have more uncertainty when you are extrapolating, the extrapolations will still be heavily influenced by the model which may or may not be very good in the extrapolated regions and there will be no way for you to tell.

Turning the model into a standard regression model

If the observations are positive and it is reasonable to consider a multiplicate rather than additive error, you can turn this problem into a standard regression problem. Then, we could assume

$$ y_{i,t} = \alpha e^{\beta/\theta_i} t e^{\epsilon_{i,t}}$$

where $e^{\epsilon_{i,t}}$ is the multiplicative error. If we take logarithms, then

$$ \log y_{i,t} = \log(\alpha) + \beta \frac{1}{\theta_i} + \epsilon_{i,t}. $$ If we assume $\epsilon_{i,t} \stackrel{ind}{\sim} N(0,\sigma^2)$, this is a standard simple linear regression model where $\log(\alpha)$ is the intercept and $\beta$ is the slope for inverse temperature. This model can be trivially fit using any regression software.

If you are using a Bayesian approach, then you can also trivially get uncertainty on $\alpha$ rather than $\log(\alpha)$ by taking samples of $\log(\alpha)$ and exponentiating. If you are using the standard prior ($p(\log(\alpha),\beta,\sigma^2) \propto 1/\sigma^2$), the forecast distribution for $\log(\tilde{y})$ is a Student $t$ distribution and can be found in most Bayesian textbooks. The statement about extrapolation above still applies.

If it is reasonable, I would certainly opt for this second approach.


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.