Interpreting a generalised linear mixed model with binomial data I have a generalised linear mixed model with binomial response data, the model:
model <- glmer(RespYN ~ Treatment + Gender + Length + (1 | Anim_ID),
  data = animDat,
  family = binomial(link = "logit"))

I am no statistician (I'm a biologist) so I have no idea how to interpret the data. With a linear mixed model I understand, due to the mean differences, etc. being means of the response variables. With binomial GLMMs I am unsure.
How do I prove that the treatment is causing/not causing the response?
Here is my output (sorry for dumping all of it):
Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
 Family: binomial  ( logit )
Formula: RespYN ~ Treatment + Gender + Length + (1 | Anim_ID)
   Data: animDat

     AIC      BIC   logLik deviance df.resid 
   142.1    158.1    -66.1    132.1      176 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.6913 -0.3021 -0.2059  0.4435  3.8066 

Random effects:
 Groups       Name        Variance Std.Dev.
 Cockroach_ID (Intercept) 0        0       
Number of obs: 181, groups:  Cockroach_ID, 10

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  5.03046    3.62723   1.387    0.165    
TreatmentPo -4.06399    0.48900  -8.311   <2e-16 ***
GenderM      0.13323    0.49365   0.270    0.787    
Length      -0.05896    0.05758  -1.024    0.306    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
            (Intr) TrtmnP GendrM
TreatmentPo -0.175              
GenderM     -0.412 -0.039       
Length      -0.995  0.139  0.348
convergence code: 0

 A: A couple of extra notes on top of what @RobertLong already answered:


*

*As Robert also noted, the interpretation of the coefficients from generalized linear mixed models are conditional on the random effects. Most often this is not the interpretation you are looking for. For more info on this check here.

*You have fitted the model with the default Laplace approximation. This may be inaccurate, especially for dichotomous data. You should better fit the model with the adaptive Gaussian quadrature by specifying a value much greater than one for the nAGQ argument of glmer().

A: The interpretation is the same as for a generalised linear model, except that the estimates of the fixed effects are conditional on the random effects.
Since this is a generalized linear mixed model, the coefficient estimates are not interpreted in the same way as for a linear model. In this case you have a binary outcome with a logit link, so the raw estimates are on the log-odds scale.
The estimated coefficient for the intercept, 5.03046, is the log odds of RespYN being 1 (or whatever non-reference value it is coded as) when Lengthis equal to zero, and  Treatment and Gender take their reference value. A value of zero for dLength might not make sense in your sample, since presumably it will never be negative and is always far above zero, and if so, you might want to consider centering it so that a zero value for the centered variable is more meaningful.
The estimate for Length of -0.05896 means that a 1 unit increase in Length is associated with a 0.05896 decrease in the log-odds of RespYN being 1, compared to RespYN being 0. If we exponentiate this number then we obtain the odds ratio of 0.9427445, which means that for a 1 unit increase in Length we expect to see (approximately) a 6% decrease in the odds of RespYN being 1.
The estimate for TreatmentPo of -4.06399 means that Treatment = Po is associated with 4.06399 lower log-odds than the other treatment group of RespYN being 1, compared to RespYN being 0. This can be exponentiated as above to obtain an odds ratio. The same analysis applies to Gender.

How do I prove that the treatment is causing/not causing the response?

Nothing can be proven with statistics, especially with observational studies. You can say that, while controlling for Gender,Length and the repeated measures within Anim_ID, you have found evidence that the association of Treatment with the the outcome is not zero. You could also say that, if the association of Treatment with the the outcome is actually zero, then the probability of observing the data that you have, or data more extreme, is less than 0.0000000000000002
Lastly, I notice that you have specified random intercepts for Anim_ID in your model formula, yet the model output says that Cockroach_ID is the grouping variable. This is rather odd, normally they would be the same. Moreover, the convergence code is zero which indicates that the model has not converged, and the estimated variance for the random effect is zero. This means that potentially there is no variation within Anim_ID. It would be a good idea to fit a model with glm() (ie without random intercepts but with Anim_ID as a fixed effect) and see how the model estimates compare.
