SemEff using gls models

[carvazgon90]

Hello!

First of all, thank you for creating this wonderful R package 😊.

I have fitted the following psem, which consists of a list of phylogenetic least squared models (gls):

psem(gls(genome_size ~ temperature + precipitation, tree.y, correlation = corBrownian(phy= newT2)),
gls(PC1_def ~ genome_size + temperature + precipitation, data.y, correlation = corBrownian(phy= newT2)),
gls(PC2_def ~ genome_size + temperature + precipitation, data.y, correlation = corBrownian(phy= newT2)),
gls(PC3_def ~ genome_size + temperature + precipitation, data.y, correlation = corBrownian(phy= newT2)),
gls(herbivory ~ PC1_def + PC2_def + PC3_def + genome_size + temperature + precipitation, data.y, correlation = corBrownian(phy= newT2)))

However, when I tried to obtain total, direct, and indirect effects using 'semEff()', I encountered the following error:

Error in getData(m) :
'data' does not contain all variables used to fit model.

I assume this error occurs because 'gls' models require both the data and a correlation structure based on a phylogenetic tree, which is not included in the 'data' object.

Is there any way to resolve this issue? I would greatly appreciate your help.

Thanks!

Carla

[murphymv]

Hi Carla,

Yes it looks like you're right. For bootstrapping, semEff requires all variables to be in 'data'. I could have a look at the code to see if I can include the correlation structure in the data prior to resampling. I'll let you know.

Thanks for highlighting the issue.

Cheers,
Mark

[murphymv]

Is the variable 'newT2' in 'data.y'? If so, I don't think there should be an issue as the function getData() does actually evaluate whatever's in the correlation argument of gls() (having looked at the code).

[carvazgon90]

Dear Mark, Thank you very much for your help. When fitting PSEMs using the gls() function, the phylogenetic tree (newT2) is not within the same data object; instead, it exists as a separate object. I also attempted to fit PSEMs using another function for phylogenetic models, the pgls() function, which takes a combined object as input, including both data and the phylogenetic tree. However, when I try to run semEff I get the following error: Error in as.data.frame.default(x[[i]], optional = TRUE, stringsAsFactors = stringsAsFactors) : cannot coerce class ‘"comparative.data"’ to a data.frame The code for psems with pgls is as follows (using piecewiseSEM v. 1.2.1) phylogeny.y <- comparative.data(newT2, data.y, SP, vcv=TRUE, vcv.dim=3) list.pgls.y = list( pgls(genome_size ~ temperature + precipitation, phylogeny.y, lambda='ML'), pgls(PC1_def ~ genome_size + temperature + precipitation, phylogeny.y, lambda='ML'), pgls(PC2_def ~ genome_size + temperature + precipitation, phylogeny.y, lambda='ML'), pgls(PC3_def ~ genome_size + temperature + precipitation, phylogeny.y, lambda='ML'), pgls(herbivory ~ PC1_def + PC2_def + PC3_def + genome_size + temperature + precipitation, phylogeny.y, lambda='ML')) psem.pgls.y = sem.fit(list.pgls.y, phylogeny.y) Thanks again for you help, respectfully, Carla Carla Vázquez-González Postdoctoral researcher (GAIN - Fulbright) twitter.com/CarlaVzquezGon1<http://twitter.com/CarlaVzquezGon1> researchgate.net/profile/Carla-Vazquez-Gonzalez<http://researchgate.net/profile/Carla-Vazquez-Gonzalez> *Department of Ecology and Evolutionary Biology - Mooney's Lab University of California Irvine. CA. USA. https://tritrophic.weebly.com/ *Evolutionary Ecology of Plant-Herbivore Interactions Misión Biológica de Galicia (CSIC) Galicia, Spain www.plantherbivory.weebly.com<http://www.plantherbivory.weebly.com/>

…

________________________________ De: Mark Murphy ***@***.***> Enviado: domingo, 1 de octubre de 2023 10:39 Para: murphymv/semEff ***@***.***> Cc: Carla Vazquez Gonzalez ***@***.***>; Author ***@***.***> Asunto: Re: [murphymv/semEff] SemEff using gls models (Issue #53) Is the variable 'newT2' in 'data.y'? If so, I don't think there should be an issue as the function getData() evaluates whatever's in the correlation argument of gls(). — Reply to this email directly, view it on GitHub<https://urldefense.com/v3/__https://github.com/murphymv/semEff/issues/53*issuecomment-1742148957__;Iw!!CzAuKJ42GuquVTTmVmPViYEvSg!NmZPZ_4Hsue9Mp7buHUbcNYaGKmtva6Re_YD9yRN0eDNJAhh0su5K_rT98X8wrbr6Sg4yAJ_eSWcjFrVLA2-TaHF$>, or unsubscribe<https://urldefense.com/v3/__https://github.com/notifications/unsubscribe-auth/BC2ENRPZBMHOAHETJI4GP5LX5GTEFANCNFSM6AAAAAA5GVRJ7A__;!!CzAuKJ42GuquVTTmVmPViYEvSg!NmZPZ_4Hsue9Mp7buHUbcNYaGKmtva6Re_YD9yRN0eDNJAhh0su5K_rT98X8wrbr6Sg4yAJ_eSWcjFrVLFb6JgmX$>. You are receiving this because you authored the thread.Message ID: ***@***.***>

[dumaskvn]

Good evening dear all,
I'm upping this issue as I face it as well, and I was wondering whether a solution had been found ?

So far I traced the
Error in getData(m) : 'data' does not contain all variables used to fit model.
error to line 105 of stdEff-fun.R file, were the "vn" variable contain the name of the phyogenetic tree used to specify the correlation, which is unlikely to be found in the model data.frame data. Removing the name of this tree from "vn" solved the issue and made the "getData" function to work.

my modifictation :
line 105 : if (!all(!vn[grepl("tree.name",vn)] %in% names(d)))

It seems to me that the other functions of this package then rely on the "update" function to compute metrics on models, which function replicate the correlation structure specified at the origin in gls models.
However, I can't get the bootstrapping function to work on gls models, and my understanding of the functions is limited..

line 352 of bootEff-fun.R
BE <- rMapply(bootEff, m, w, SIMPLIFY = FALSE)

rMapply(bootEff, m, w, SIMPLIFY = FALSE)

ORDINARY NONPARAMETRIC BOOTSTRAP

Call:
boot::boot(data = x, statistic = s, R = R, parallel = parallel,
ncpus = nc, cl = cl)

Bootstrap Statistics :
WARNING: All values of t1* are NA
Warning messages:
1: In (function (m, w) :
101 model fit(s) or parameter estimation(s) failed. NAs reported/generated.
2: In attributes(B$t)[c("sim", "seed", "n")] <- c(B$sim, seed, nrow(B$data)) :
number of items to replace is not a multiple of replacement length

Could you explain what is going on here ?

Sincerely yours,
Keyvan

[murphymv]

Hi,

These issues with GLS models should be fixed in the latest version of the package (0.7.0 - on CRAN soon but can test with dev version in meantime). If you get a chance to test things again let me know how it goes (and sorry delays in responding).

Cheers,
Mark

[dumaskvn]

Hi !
Thanks for the reply ;)
For what I'm concerned I worked around the issue by not using GLS models (which did a worse fit than classical glm anyway in my case), so I'm I will not be able to give you feedbacks on this.
Thanks again for SemEff, it does such a good work !

[murphymv]

No prob, figured you'd find a way around it somehow! The fixes were fairly simple in hindsight and had to do with the evaluation of the correlation structure objects etc. that you supply to the GLS models. Anyway thanks again for using the package, I'm just sorry I can't maintain it as regularly as previously.

Cheers,
Mark