Korrelierte Kovariablen und Gamma Verteilung

[FledermausR]

Hallo,
ich untersuche den Einfluss von menschlicher Anwesenheit auf die Flugaktivität von Fledermäusen im Winterquartier. Ich habe folgende Variablen

Treatment (=zweifaktorieller fixed Effect (Mensch anwesend / Mensch nicht Anwesend))
Flugaktivität (FA) (= response variable, Zählung, wie oft eine FM durch die Lichtschranke geflogen ist)
Bunkertemperatur (Tb) (= numerische covariable, da FM bei höheren Temperaturen generell aktiver sind)
Außentemperatur (To) (= " " )
Bunker ID(= Random effect)

In den wärmeren Bunkern sitzen aber weniger FM (1-6) als in den kälteren Bunkern (40-90). Deshalb habe ich die Anzahl der Flugbewegungen durch die Anzahl der FM in den entsprechenden Bunkern geteilt, um das Verhältnis der Aktivität im Bezug auf die Bunkertemperatur korrekt wiederzuspiegeln. Jetzt habe ich keine ganzen Zahlen mehr und kann die Poisson Verteilung bei meinem glmer nicht mehr anwenden. Zwei Lösungsansätze, zu denen ich eine Meinung bräuchte:

1. Da To und Tb korreliert sind, würde ich Tb am liebsten weg lassen. Dann könnte ich die Flugaktivität als response variable nehmen, ohne sie durch die Anzahl der Individuen zu teilen, hätte ganze Zahlen und könnte die Poisson Verteilung nehmen. Der Korrlationskoeffizient liegt aber nur bei 0.25, und das ist glaube ich nicht groß genug, um den Rausschmiss von Tb zu rechtfertigen. Wobei man sagen muss, dass die Bunkertemperaturen relativ konstant sind und nur leicht schwanken im Vergleich zu den Außentemperaturen und die Korrelation deshalb evtl so schwach ist?!

2. Ich habe versucht ein glmer mit einer Gamma Verteilung zu erstellen.

Code: Alles auswählen

m1<-glmer(FA.count.r~treatment +  To + Tb + (treatment:To)+ (treatment: Tb) +  (1|bunker) , data=model, family="Gamma"(link="log"))

bekomme dann folgende Fehlermeldung

Code: Alles auswählen

Error in eval(family$initialize, rho) : 
  non-positive values not allowed for the 'gamma' family

Habe aber keine negativen Werte in meiner response variable!

Ich bin dankbar über jeden Hinweis.
Vielen Dank

[Athomas]

non-positive values not allowed for the 'gamma' family

Habe aber keine negativen Werte in meiner response variable!

Hmmm, hmmm - was könnte das denn sein, weder positiv noch negativ?
Null Ideen !

[EDi]

Ich stimmt Athomas zu...

Du können ein "zero-inflated Gamma" (ZIG) versuchen (z.b. mit gamlss, glmmTMB oder brms) - das ist aber nochmal eine Schippe komplexer...

In den wärmeren Bunkern sitzen aber weniger FM (1-6) als in den kälteren Bunkern (40-90). Deshalb habe ich die Anzahl der Flugbewegungen durch die Anzahl der FM in den entsprechenden Bunkern geteilt, um das Verhältnis der Aktivität im Bezug auf die Bunkertemperatur korrekt wiederzuspiegeln. Jetzt habe ich keine ganzen Zahlen mehr und kann die Poisson Verteilung bei meinem glmer nicht mehr anwenden. Zwei Lösungsansätze, zu denen ich eine Meinung bräuchte:

Hier kann ich nicht ganz folgen....
Du hast doch den random-intercept im Model - das macht der doch schon (also unterschiede zwischen Bunkern aufgreifen).

Ich würde die ganzzahlen lassen, den bunker random intercept drin lassen und eine Negativ-binomial Verteilung (oder quasi-poisson oder log-transform) ansetzen (Poisson passt nie mit ökologischen Daten - habs zumindest bisher nicht gesehen).


Mit den zwei Temperaturen müsste man nochmal nachdenken...

Hast du mehrmals an den Bunkern gemessen? Dann könnte man eine Temperatur noch mit in den Random-part reinnehmen um den Bunker effekt feiner aufzuschlüsseln (z.b. also random slope pro bunker - da brauchst du aber einiges an Messungen pro bunker!) - macht das Sinn? (ist schon spät und ohne das Design zu kennen ist das rumwunderei....).

Du könntest auch jeweils nur eine reinnehmen - die zwei Modelle vergleichen und das bessere nehmen.
Oder du sagts, dein Bunker random-intercept greift Tb sowieso schon auf und lässt das a-priori raus.

[FledermausR]

Vielen Dank

Ich würde das Model generell gerne so einfach wie möglich halten, deshalb würde ich gerne versuchen, ohne ein "zero inflated Gamma" Model auszukommen

Zum Versuchsaufbau:
Ich war neun mal in jedem der acht Bunker (in einem wöchentlichen Abstand). Als Kontrollwert für einen Tag ohne menschlichen Besuch habe ich die Aktivität von dem jeweiligen Bunker am Tag vor meinem Besuch gewählt. Ich habe also für jeden Bunker je 9 Datenpunkte für das Treatment "monitoring" und "control". Für jeden Tag (sowohl für den Kontrolltag als auch für den Monitoringtag) habe ich für jeden Bunker eine Messung von der Innentemperatur des Bunkers.

Hast du mehrmals an den Bunkern gemessen? Dann könnte man eine Temperatur noch mit in den Random-part reinnehmen um den Bunker effekt feiner aufzuschlüsseln (z.b. also random slope pro bunker - da brauchst du aber einiges an Messungen pro bunker!) - macht das Sinn? (ist schon spät und ohne das Design zu kennen ist das rumwunderei....).

Finde ich gut! Habe ich ausprobiert:

Code: Alles auswählen

m1<-glmer(FA~treatment * To + (1+Tb|bunker), data=model, family="poisson")

Kann ich das so machen, oder muss Tb als main effect/Kovariable mit in das Model aufgenommen werden? Geht das überhaupt so? Weil meines Wissens nach ein random effect eine faktorielle Variable sein muss, was Tb aber nicht ist.

der output wäre folgender:

Code: Alles auswählen

> Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: poisson  ( log )
Formula: FA ~ treatment * To + (1 + Tb | bunker)
   Data: model

     AIC      BIC   logLik deviance df.resid 
  1335.7   1356.2   -660.8   1321.7      132 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.7968 -1.7405 -0.3924  1.1938 11.8905 

Random effects:
 Groups Name        Variance Std.Dev. Corr 
 bunker (Intercept) 1.08965  1.0439        
        Tb          0.02358  0.1536   -0.97
Number of obs: 139, groups:  bunker, 8

Fixed effects:
                 Estimate Std. Error z value Pr(>|z|)    
(Intercept)      1.844755   0.105661  17.459  < 2e-16 ***
treatmentmon     0.461424   0.058851   7.841 4.49e-15 ***
To               0.104182   0.009664  10.781  < 2e-16 ***
treatmentmon:To -0.069205   0.012463  -5.553 2.81e-08 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
            (Intr) trtmnt To    
treatmentmn -0.342              
To          -0.234  0.374       
tretmntmn:T  0.169 -0.361 -0.765

Vielen Dank !!

[EDi]

Geht das überhaupt so? Weil meines Wissens nach ein random effect eine faktorielle Variable sein muss, was Tb aber nicht ist.

Ist dann halt ein Random slope & random intercept model.