[R-sig-ME] R-sig-mixed-models Digest, Vol 175, Issue 23

Tue Jul 27 19:54:30 CEST 2021

Hello,
I have a question which will probably sound very basic to most of you. However, I am insecure about my results and rather ask before presenting and publishing something I am not 100% certain about. Therefore, I would really appreciate your help!

I am working on an ecotoxicological field study where a fungicide was applied onto 4 fields (with 7 subplots) whereas nothing was applied onto 4 control fields. The experiment takes place over 43 days. I checked for spatial and temporal autocorrelation.  My dependent variable is for this example the number of insects. I am unfortunately not allowed to provide data, but I believe that my questions can be answered nevertheless.
DAT is a continuous variable (days of the experiment)
Treatment is either TRUE or FALSE
I use a linear model, because the number of insects is over 2000

I constructed this model: MyModel <-  glmmTMB(NumberInsects ~ Treatment poly(DAT,2) + (1|Plot/Subplot), data=dat, REML=T, family=gaussian(link="identity"))

I checked for the model fit via Dharma residual diagnostics, including temporal and spatial correlations.

I answered the question for the significance of difference of number of insects in respect to treatment and control on specific days via:
A <- emmeans(MyModel, ~Treatment|DAT, at =list(DAT=c(0,8,14,28,35,43)))
pairs(A)

However, I have issues answering the following question:

  *   Does treatment have an effect on the number of insects, considering the complete time course of the experiment?

Is it correct to use the emtrends approach here? I understand that the slope of the specified time frame are compared, which would also decouple the results from the initial numer of insects, as the slope is independent from the intercept.
This is the code that I am using to answer this question:
B <- emtrends(MyModel, "Treatment", var="DAT")
pairs(B)

It would be wrong to answer this question with the p-value of Treatment of the summary(myModel) function, right? However, I am wondering whether the 2nd polynomial of DAT is regarded in the emtrends approach.

I really hope that you can help me!

Kind regards,

InkaMarei

________________________________
Von: R-sig-mixed-models <r-sig-mixed-models-bounces using r-project.org> im Auftrag von r-sig-mixed-models-request using r-project.org <r-sig-mixed-models-request using r-project.org>
Gesendet: Dienstag, 27. Juli 2021 12:00
An: r-sig-mixed-models using r-project.org <r-sig-mixed-models using r-project.org>
Betreff: R-sig-mixed-models Digest, Vol 175, Issue 22

Send R-sig-mixed-models mailing list submissions to
        r-sig-mixed-models using r-project.org

To subscribe or unsubscribe via the World Wide Web, visit
        https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
or, via email, send a message with subject or body 'help' to
        r-sig-mixed-models-request using r-project.org

You can reach the person managing the list at
        r-sig-mixed-models-owner using r-project.org

When replying, please edit your Subject line so it is more specific
than "Re: Contents of R-sig-mixed-models digest..."

Today's Topics:

   1. Re: Multilevel equation (Brian Hudson)

----------------------------------------------------------------------

Message: 1
Date: Mon, 26 Jul 2021 22:05:19 -0400
From: Brian Hudson <bhudson.gsu using gmail.com>
To: Thierry Onkelinx <thierry.onkelinx using inbo.be>
Cc: r-sig-mixed-models <r-sig-mixed-models using r-project.org>
Subject: Re: [R-sig-ME] Multilevel equation
Message-ID:
        <CANodPHM1Rg4c-B7x_whDewOcy6mm_z09T83zTmxom=yRSZozSA using mail.gmail.com>
Content-Type: text/plain; charset="utf-8"

Thierry,

Thank you! I appreciate your help and explanation - that makes sense and I
can see where my other attempts were incorrect.

A couple questions-

1) The two lines that defined eta confused me - are they supposed to be
equal to each other? I edited the equation below such that pi has the link
function and eta has the linear equation - does that work?
2) I added epsilon for the individual error
3) There was no k index (individual level) in the equations, just i and j,
so i added some indexing in the predictors.

Does this equation make sense? Any issues with what I did? Thanks again for
your (and the community's) help.

https://quicklatex.com/cache3/9f/ql_4a4eb44285f65ea0dcaf93d551c44c9f_l3.png

$$b_i\sim \mathcal{N}(0, \sigma_s^2)$$
$$b_{ij}\sim \mathcal{N}(0, \sigma_{y}^2)$$
$$\eta_{ijk} = \beta_0 + \beta_1 \textrm{X}\textsubscript{m[i]} + \beta_2
\textrm{X}\textsubscript{y[i,j]} + \beta_3
\textrm{X}\textsubscript{s[i,j,k]} + b_i + b_{ij} + \epsilon_{ijk}$$
$$\pi_{ijk} =\frac{e_{ijk}^{\eta}}{1+e_{ijk}^{\eta}}$$
$$Y_{ijk} \sim Binom(1, \pi_{ijk})$$

On Mon, Jul 19, 2021 at 1:34 PM Thierry Onkelinx <thierry.onkelinx using inbo.be>
wrote:

> Dear Brian,
>
> I'd write it as follows. In the case of a Gaussian model, you only have to
> write $Y_{ijk} \sim \mathcal{N}(\eta_{ijk}, \sigma^2)$ and drop the link
> function. (And you could replace \eta with \mu). Basically, Y depends on a
> distribution defined by some parameters. And these parameters might need
> some further definition.
>
> $i$: state index
> $j$: year index
> $k$: observation index
> $X_m$: state_mnthyr_pred
> $X_y$: state_year_pred
> $X_s$: state_pred
> $$Y_{ijk} \sim Binom(\pi_{ijk})$$
> $$\eta_{ijk} = \frac{\pi_{ijk}}{1- \pi_{ijk}}$$
> $$\eta_{ijk} = \beta_0 + \beta_1X_m + \beta_2 X_y + \beta_3 X_s + b_i +
> b_{ij}$$
> $$b_i\sim \mathcal{N}(0, \sigma_s^2)$$
> $$b_{ij}\sim \mathcal{N}(0, \sigma_{y}^2)$$
>
> Best regards,
>
>
> ir. Thierry Onkelinx
> Statisticus / Statistician
>
> Vlaamse Overheid / Government of Flanders
> INSTITUUT VOOR NATUUR- EN BOSONDERZOEK / RESEARCH INSTITUTE FOR NATURE AND
> FOREST
> Team Biometrie & Kwaliteitszorg / Team Biometrics & Quality Assurance
> thierry.onkelinx using inbo.be
> Havenlaan 88 bus 73, 1000 Brussel
> www.inbo.be<http://www.inbo.be>
>
>
> ///////////////////////////////////////////////////////////////////////////////////////////
> To call in the statistician after the experiment is done may be no more
> than asking him to perform a post-mortem examination: he may be able to say
> what the experiment died of. ~ Sir Ronald Aylmer Fisher
> The plural of anecdote is not data. ~ Roger Brinner
> The combination of some data and an aching desire for an answer does not
> ensure that a reasonable answer can be extracted from a given body of data.
> ~ John Tukey
>
> ///////////////////////////////////////////////////////////////////////////////////////////
>
> <https://www.inbo.be>
>
>
> Op ma 19 jul. 2021 om 17:44 schreef Brian Hudson <bhudson.gsu using gmail.com>:
>
>> Hello,
>>
>> I am fitting a multilevel model in `lme4` and am having trouble writing
>> the
>> equation for it. I very much appreciate any help. The formula and code is
>> below, but I am not sure if the equation represents the error correctly -
>> do I need to include error terms or is that captured by the distributions?
>> I am also not sure if I am representing the logit function correctly with
>> the indexing or functional form.
>>
>> The data are comprised of US-state months nested within US-state-years and
>> US-states. I include predictors at each level and a varying intercept for
>> both state-years and states.
>>
>> The formula looks like this in R:
>>
>> ```
>> as.formula(outcome ~ state_mnthyr_pred + state_year_pred + state_pred +
>>                          (1 | state) + (1 | state_year))
>> ```
>> Where the outcome is dichotomous. The state months (e.g. jan-2010,
>> feb-2010
>> ... jan-2013) are nested with state years and within states.
>>
>> The formula I am using can be seen here:
>>
>>
>> https://quicklatex.com/cache3/e9/ql_038eeb4e4e1b0af94d3ef69fe4ff7be9_l3.png
>> And the LaTeX code:
>>
>> $$
>> \begin{aligned}
>>     \mu &=\alpha_{j[i],k[i]} +
>> \beta_{0}(\operatorname{state\_mnthyr\_pred})\ \\
>>     \alpha_{j}  &\sim N \left(\gamma_{0}^{\alpha} +
>> \gamma_{1}^{\alpha}(\operatorname{\textrm{state\_year\_pred}}),
>> \sigma^2_{\alpha_{j}} \right)
>>     \text{, for \textrm{State-Year} j = 1,} \dots \text{, J} \\
>>     \alpha_{k}  &\sim N \left(\gamma_{0}^{\alpha} +
>> \gamma_{1}^{\alpha}(\operatorname{\textrm{state\_pred}}),
>> \sigma^2_{\alpha_{k}} \right)
>>     \text{, for State k = 1,} \dots \text{, K}\\
>> \pi_{i} &=\frac{e_{i}^{\mu}}{1+e_{i}^{\mu}}\\
>> y_{i j k} \sim & \operatorname{Binom}\left(1, \pi_{i}\right)\\
>> \end{aligned}
>> $$
>>
>> I really appreciate any help. Thank you.
>>
>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-mixed-models using r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>

        [[alternative HTML version deleted]]

------------------------------

Subject: Digest Footer

_______________________________________________
R-sig-mixed-models mailing list
R-sig-mixed-models using r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

------------------------------

End of R-sig-mixed-models Digest, Vol 175, Issue 22
***************************************************

	[[alternative HTML version deleted]]