Time Series of Wound Healing percentages / proportions

Question

I am looking at wound healing over time and have hit a wall with analysis. My data consist of injury cases which are tracked over time. Injury healing is measured as percentage wound area/perimeter change since initial sighting. Each case varies in the number of days since initial sighting and the total number of sightings, with some only seen twice and some seen 5 times.

DV = Area Change

IV = Days, Severity, Location, Type

Random effect = Shark ID, Injury

I want to look at changes over time and be able to statistically explore whether wounds are decreasing in size over time and what factors can influence this. Key points:

My response is a percentage that is dependent on previous records. How can I combat this?

Here is an example of data:

Injury  Shark ID   Name      Date     Days    Severity  Location  Type        Cause      Area Change     Perimeter Change  
1       WS198      Naococco  24.6.15     0    Major     Flank     Laceration  Boat                 0                    0
1       WS198      Naococco  8.7.15     14    Major     Flank     Laceration  Boat        74.6858632           63.6698483  
1       WS198      Naococco  10.7.15    16    Major     Flank     Laceration  Boat        81.4165492           77.0354486  
1       WS198      Naococco  1.8.15     38    Major     Flank     Laceration  Boat        81.9152998            89.704627  
1       WS198      Naococco  4.8.15     41    Major     Flank     Laceration  Boat        86.2165258           89.1497158  
1       WS198      Naococco  18.8.15    55    Major     Flank     Laceration  Boat        86.1513557           90.3776043  
1       WS198      Naococco  19.9.15    87    Major     Flank     Laceration  Boat        97.9639636           93.7961375 
1       WS198      Naococco  24.9.15    92    Major     Flank     Laceration  Boat        96.8353572            91.813769     
1       WS198      Naococco  3.11.15   132    Major     Flank     Laceration  Boat        99.7190944           98.0207677 
1       WS198      Naococco  10.11.15  139    Major     Flank     Laceration  Boat        99.8686985           98.8372268  
1       WS198      Naococco  30.11.15  159    Major     Flank     Laceration  Boat               100                  100 
2       WS198      Naococco  27.3.16     0    Major     Flank     Laceration  Boat                 0                    0  
2       WS198      Naococco  29.3.16     2    Major     Flank     Laceration  Boat         28.902414           17.0816156   
2       WS198      Naococco  4.4.16      8    Major     Flank     Laceration  Boat        47.1829361           45.5672997  
2       WS198      Naococco  30.5.16    64    Major     Flank     Laceration  Boat        97.4773761           91.7192607 
2       WS198      Naococco  13.9.16   170    Major     Flank     Laceration  Boat               100                  100 
3       WS198      Naococco  5.10.17     0    Major     Flank     Laceration  Boat                 0                    0 
3       WS198      Naococco  26.10.17   21    Major     Flank     Laceration  Boat        65.2492641           67.0574269  
3       WS198      Naococco  30.10.17   25    Major     Flank     Laceration  Boat        68.9993994           66.1401071 
3       WS198      Naococco  10.12.17   66    Major     Flank     Laceration  Boat               100                  100  

Please tell me if this does not make sense or if more information is required to clarify my question.

I know I cannot use a standard lmer model, as follows, because of my response:

healing.lmer <- lmer(Area Change ~ Days + Severity * Type * Location +
                (1 | Shark ID) + (1 |Injury), data = healing)

How can I change my response or my model to fit?

Answer 1

This is a complex (and interesting) dataset and I think the reason why you haven't gotten any specific answers yet is because people would need to play around with the actual data to better understand it. I would say there is no ONE answer to your question.

I would go ahead and do some initial data visualizations as suggested by @JamesPhillips. That should give you a very good idea what's going on in the data (you most likely know this already). Then I would start fitting some models that take into account the dependency (repeated measures) in your data. You can simply start with a linear-mixed effect model as specified and do detailed model validations. I agree with you that given your outcome variable a linear-mixed effects model may not be ideal but it could be a start. What I would probably do next then using the glmmTMB package, which handles Beta mixed models including zero and zero-one inflated models. If you don't want to go down that rabbit hole, you could try a logit or arcsine transformation of your percentages, however see: The arcsine is asinine: the analysis of proportions in ecology.

Regarding your random effects: Is the injury coding a running number across your sharks or does it always start from $1,...,n$ within each shark (cannot see this in your example)? If it is the latter, then you need to explicitly nest your random effect as + (1|Shark ID / Injury) (which is the same as + (1|Shark ID) + (1|Shark ID : Injury)). If it is the former (unique running number across sharks), you can either do + (1|Shark ID / Injury) or + (1|Shark ID) + (1|Injury) as they are equivalent in this case.

Also have a look at @BenBolker 's GLMM FAQ. It pretty much answers all technical questions regarding mixed models that you could have.

Answer 2

I agree this one is challenging, but you can design dynamic models that will respect the proportional nature and the temporal dependence of your observations. The {mvgam} package in R was designed specifically for this kind of challenging dataset. The package works by formulating State-Space models where your data (usually in the form of time series) are considered be be noisy observations of some true latent dynamic process. The process is modeled using appropriate dynamic models (in your case, a continuous time autoregressive model might be appropriate given the irregular spacing of your observations in time). You can also include any necessary predictors in the process model that you think impact its evolution, and these predictors can be in the form of penalized smooth functions like you would use in Generalized Additive Models. The observations are then taken using a separate model, and this is linked to the process model using an appropriate link function. This allows you to respect the features of your observations in the same way you would by using a Generalized Linear Model. A relevant resource to showcase how you could handle proportional time series would be the {mvgam} vignette on time-varying effects.