florianhartig
commented
4 years ago Hi, I just played a bit around with that - here an example of how one could modify the residuals function in piecewiseSEM to support DHARMa ... maybe one should check in the loop if the fitted models are supported by DHARMa, and fall back to standard otherwise.
Add-on: obviously, the model doesn't make a lot of sense, and also DHARMa doesn't add much for lms, but just as a proof of principle.
library(piecewiseSEM)
mod = psem(
lm(rich ~ distance + elev + abiotic + age + hetero + firesev + cover, data = keeley),
lm(firesev ~ elev + age + cover, data = keeley),
lm(cover ~ age + elev + abiotic, data = keeley)
)
summary(mod)
plot(mod)
library(DHARMa)
residuals2 <- function(fittedPSEM, type = "DHARMa"){
if(type == "standard"){
return(residuals(fittedPSEM))
}
if(type == "DHARMa"){
sim = list()
res = list()
for(i in 1:(length(fittedPSEM) - 1)){
sim[[i]] <- simulateResiduals(fittedPSEM[[i]])
res[[i]] = sim[[i]]$scaledResiduals
}
res = matrix(unlist(res),ncol=length(res),byrow=F)
return(list(residuals = res, DHARMaSimulations = sim))
}
}
res = residuals2(mod)
res$residuals
# plots and tests for model 1
plot(res$DHARMaSimulations[[1]])
testResiduals(res$DHARMaSimulations[[1]])
# in principle, it would be possible to plot residuals also against a particular path, or test for homogeneity along this path
jslefche
jebyrnes
Hi Jon,
I'm the developer of DHARMa and incidentally just attending a workshop on SEMs, which made me think about the possibility to integrate DHARMa with piecewiseSEM.
I have looked at the structure of your package, and I see two routes one could go
Integrating DHARMa into piecewiseSEM: most of the models that are supported by piecewiseSEM could also be checked with DHARMa (except nlme). So, it would be possible to add switch in piecewiseSEM and return DHARMa residuals for each model, or p-values for DHARMa omnibus residual tests that could be superimposed on the DAG or displayed in piecewiseSEM.
Integrating piecewiseSEM into DHARMa: to support piecewiseSEM, minimum requirements are described here https://github.com/florianhartig/DHARMa/wiki/Adding-new-R-packages-to-DHARMA
Option 2 would allow to create omnibus for a tested SEM from DHARMa, e.g. global dispersion of spatial autocorrelation tests. Given the nature of the package, however, it seems to me more logical to consider each link individually (because averaging over them might mask problems), which would suggest to me that it would be more logical to go via route 1. I think what could added easily is implement
a) the ability to return DHARMa instead of normal residuals per link b) the option to impose the results of testResiduals (3 general tests) or testUniformity (tests just for distribution) on the links in the dag, e.g. in the summary() or in the plot() function, to alert users to connections that are problematic