TL;DR: Is there some way I can do 'random knots' (where the position of the knot varies according to some group, not the number)?

Full question: After obtaining experimental data, I find that one of the covariates has an interesting relationship with the response in most of the subjects: Beyond a certain threshold, all or almost all of the response are horizontal (at the maximum value). (Example 1, Example 2) I figured it would make sense to fit the final model to incorporate this visual finding. However, the position of the knot is obviously different for each subject.

Is there an existing method to do this? (Aside from this covariate and the random effect of the subject, there's also one more fixed effect and one more random effect.)

P.S. I have seen the answer to this question, but I think I have good reason to prefer piecewise linear to cubic splines here: The result will be easier to interpret (this is important as our focus is inference, not prediction) and most of the graphs look like piecewise linear will capture the relationship well enough without introducing curves.

P.P.S. One thing I have in mind is to do piecewise linear regression for each subject, determining knot position by just a grid search, then fit the residuals on everything else. I'm not sure if this is the best way to go about it though.

Edit: Currently I'm trying to treat the response as censored instead.

  • $\begingroup$ please provide more detail about how the data arise. What is the experimental design ? $\endgroup$ Jun 18 '18 at 16:04
  • $\begingroup$ @RobertLong The response variable is the judgement from human participants (from 0 to 100, normalised to be 0 to 4 in the graph for easier interpretation). The covariates (including the one on the x-axis of the graph) are two theoretical values that affect the human judgement (one of our main concerns is to find how the two are weighted, as there is no previous literature on this). The random effects are item and subject. $\endgroup$ Jun 20 '18 at 15:51

You can do this in the R package mcp. Although your actual full model may be outside the scope of mcp, this is a way to do "random effects" change points.

The mcp package contains a demo dataset called ex_varying:

> library(mcp)
> head(ex_varying)
    id  x id_numeric         y
1 John  1          5 30.792018
2 John  5          5  1.027091
3 John  9          5 58.793870
4 John 13          5 40.300737
5 John 17          5 57.566408
6 John 21          5 80.876520

Model two joined slopes with the change point location varying by id. You will recognize this syntax from lme4:

model = list(
  y ~ 1 + x,          # intercept + slope
  1 + (1|id) ~ 0 + x  # joined slope, varying by id
fit = mcp(model, ex_varying)
plot(fit, facet_by = "id", cp_dens = FALSE)

enter image description here

You can visualize the change point posteriors using plot_pars(fit, "varying") and summarise them using ranef(fit). Read more in the mcp article on random effects (called "varying effects" in mcp cf. the terminology from the brms package).

  • $\begingroup$ Thanks a lot! The project I was asking the question for is already done and I did the analysis in a rather different way, but I'll definitely keep mcp in mind if something else comes up in the future. $\endgroup$ Jan 19 '20 at 20:18
  • $\begingroup$ @WavesWashSands Great! It is good practice to accept answers that solve the problem you posed to help future users find relevant answers - no matter whether your "real" problem changed or you found some other solution. If you found a better solution yourself, you can also post that and accept your own answer. $\endgroup$ Jan 20 '20 at 9:55

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.