Representing uncertainty along a curve (error bands)

Interfaces brms
1
  • Operating System: Windows10
  • brms Version:2.8

I created a nonlinear mixed-effects model.
Now I would like to generate a plot with my actual data points and a curve with parameter values, which are the parameter estimates obtained when fitting the model.
However, if I only plot these, I lose a lot of information about the uncertainty of the parameters.
So I was wondering what the best way is to represent this uncertainty?

A) Should I just draw samples from the distributions of the parameters (they all have normal prior, so maybe using multivariate normal?) and, generate, say, 100 curves with these, and plot them at the same plot as my fitted curve?

B) Or is there any way to plot error bands to represent that uncertainty using brms?
Or should I just use the estimate I got for Sigma (in my case family=gaussian()) ´,
even though in this case I would lose information about the uncertainty of Sigma?

I would really appreciate it if someone could give me some insight.

2

I am not sure I understand what you want to do, because your title mentions “fitted curve” whereas you speak of parameters in the text.

Do you want a plot with one predictor on the x-axis, and the criterion/outcome on the y-axis? Then the marginal_effects function might be of use. You could use marginal_effects to plot predictions with confidence bands, and then manually add the observed data as points.

More generally about uncertainty:
As far as I can tell, confidence bands typically show confidence/credible intervals for the mean, without taking into account Sigma. If one takes Sigma into account, one often speaks of “prediction intervals”. marginal_effects internally uses the posterior_predict functions, and thus takes Sigma (including uncertainty about Sigma) into account.

3

Hi, Guido,

Thank you for your reply! I tried to edit the title so it’s a bit more descriptive and corresponding to the content of my post.
I looked into the marginal_effects function, it seems useful.
But in this case it seems to be generating a plot of a curve which takes into account all data.
What I would like to do is generate a marginal_effects plot for each subject (“Sample”) in my dataset (each subject has several measurements). I am assuming this somehow takes into account the value of the random effects for each specific subject. However, I am not quite sure how to do it.
I tried

marginal_effects(fit, conditions=dataset$Sample, re_formula = NULL, method=“predict”)

But this only resulted in R plotting M identical plots, where M is the number of all DATA POINTS in my dataset, rather than the number of distinct subjects (in my case named “Sample”).

I also tried

newdata <-data[ which(data$Sample==1),]
marginal_effects(fit, conditions=newdata, re_formula = NULL, method=“predict” ),

Which resulted in N identical plots, where N is the number of data points, corresponding to subject #1.

Have you had experience with a problem similar to mine?
If so, how did you manage to plot a separate marginal_effects plot for each subject?

4

I also want to obtain such upper and lower band curves for ROC like curves.

My poor idea is the followings, but I am not sure whether it is correct or not. I also want to know how to draw such band curves.

In case that it is impossible to calculate by formula
If model is analytically difficult to calculate, then sampling is the only way.
For example, by bootstrap sampling or jacknife sampling, we can replicate various 100 datasets from the empirical distribution.

Then, we also obtain 100 curves,
\gamma^1, \gamma^2,...,\gamma^{100}, \gamma^i (t)=(t,y^i(t)) \in \mathbb{R}, i=1,2,...,100 for each replicated dataset.

Then for each t, we can evaluate the “95% interval” of the data \{ y^1(t),y^2(t),...,y^{100}(t) \} under approximate assumption that the data \{ y^1(t),y^2(t),...,y^{100}(t) \} are normally distributed for each t.

In case that it is possible to calculate by formula

The most traditional way is the following.
Let f(y|\theta) be a likelihood of model. Then curve can be written using model parameter \theta. \gamma_\theta :\mathbb{R} \to \mathbb{R}^2 ;t \mapsto (t,y _{\theta}(t)). If a dataset D is given then using poterior mean \hat{\theta}(D) we get the curve (t,y _{\hat{\theta}(D)}(t)) for each t.

Using posterior predictive distribution, which is measure on the set of all datasets \{D\}, we can calculate the variance of the curve (t,y _{\hat{\theta}(D)}(t)) for each t and it is the desired band curves.

5

I am still not 100% sure what you want to do. In your title it reads “parameter estimates” and in the text it sounds as if you want to plot predictions from the model (which is also consistent with your use of the predict function). In my understanding, a parameter is a parameter of the model, which only together with the data makes a prediction.

Anyhow, I am still assuming you want to plot predictions. In particular, you seem to want one plot per participant with a predictor/conditions on the x-axis and observed/predicted values on the y axis*. Correct?

I do not know if brms has a function to do this automatically. You could try the ggeffects package, but I am not sure this has a function that makes it easy to make the plot you (seem to) want to make.

To make the plot “manually”, you can use the posterior_predict function in brms to generate predictions. This will return an number_of_sampels * number_of_observations matrix (hope I am getting the dimensions right here…) with posterior predictions, from which you can calculate mean/median and credible intervals for any group of observations you wish.
Once you have these data, you can use base R or ggplot to put the plot you want together.

*it is not clear if the predictor that varies is nominal/ordered/continuous. Depending on what it is, you might want to plot points with CIs (nominal, ordered) or lines with confidence bands (continuous).

6

Maybe you want to have a look at the tidybayes package: Extracting and visualizing tidy draws from brms models • tidybayes
and especially the spread_draws() and the add_predicted_draws()/add_fitted_draws() functions.

Regards!

3 Likes