Beta estimates outside of -3 to 3 I am performing logistic regression with glmmTMB with random intercepts and random slopes in R, with 5 fixed effects. I have 9 different datasets that I am running through the same model. Before running my model I have centered and scaled each dataset. If I understand correctly, my beta estimates should be between -3 and 3 because of the standard deviations associated with scaled and centered data. Well, in 4 of my datasets, the beta estimates for one particular covariate range from -3.13 to -5.89. I have double checked these datasets and everything appears to be fine. Depending on the dataset, my sample size is 28-32 animals.
So my questions are:

*

*Is my understanding correct that beta estimates should be between -3 and 3 with scaled and centered data?

*If this is true, what could be causing my beta estimate to be outside of this range?

*Could it be due to too much (or too little) variation for this particular covariate? If so, what can I do about this?

 A: I'll give an example using a vanilla logistic regression (no random effects) to show why this assertion that the coefficients must between $-3$ and $3$ need not hold.
set.seed(2022)
N <- 1000
x1 <- rnorm(N, 0, 1)
x2 <- rnorm(N, 0, 1)
z <- 4*x1 - 4*x2 # Linear combination of the x1 and x2 features
pr <- 1/(1 + exp(-z)) # Transform through the inverse logit
y <- rbinom(N, 1, pr) # Draw from the conditional binomial distributions
x1s <- (x1 - mean(x1))/sd(x1) # Standardize x1
x2s <- (x2 - mean(x2))/sd(x2) # Standardize x2
Ls <- glm(y ~ x1s + x2s, family = binomial) # Fit with standardized features 
summary(Ls)

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)   0.1534     0.1230   1.246    0.213    
x1s           3.7805     0.2878  13.138   <2e-16 ***
x2s          -4.3728     0.3215 -13.599   <2e-16 ***

Interestingly, the coefficients are not that different when we fit to the original features.
L <- glm( y ~ x1  + x2 , family = binomial) # Fit with original features
summary(L)

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)   0.1328     0.1228   1.081     0.28    
x1            3.7872     0.2883  13.138   <2e-16 ***
x2           -4.2710     0.3141 -13.599   <2e-16 ***

A: Linear regression or logistic regression are linear models, scaling is a linear transformation. Generalized linear models take the form
$$\begin{align}
\eta &= \beta_0 + \beta_1 X_1 + \dots + \beta_k X_k \\
E[y|\mathbf{X}] &= g(\eta)
\end{align}$$
After scaling, we want to get the same, optimal, result as before, so we want $g(\eta)$ to be the same, so $\eta$ needs to be the same. That means, we can focus only on the linear predictor $\eta$ and the reasoning would be the same regardless of the link function and family of the GLM.
We use scaling of the data for two reasons: to improve interpretability of the parameters, and in some models (e.g. regularized) to prevent numerical problems in optimization. Notice that when you scale the data, the result would be the same as with non-scaled data, but the parameters would adapt to the scaling. So if you have linear function $\alpha + \beta X$ then after scaling $X' = (X-\bar x)/s_X$ you have to also scale the parameters accordingly. So you have
$$
\require{cancel} 
\begin{align}
\alpha + \beta X &= \alpha' + \beta' X'\\
&= (\alpha + \beta' \tfrac{\bar X}{s_X}) + \beta' \tfrac{X - \bar X}{s_X} \\
&= (\alpha + (\beta s_X) \tfrac{\bar X}{s_X}) + (\beta s_X) \tfrac{X - \bar X}{s_X} \\
&= \alpha + \beta \cancel{s_X} \tfrac{\bar X}{\cancel{s_X}} + \beta \cancel{s_X} \tfrac{X}{\cancel{s_X}} - \beta \cancel{s_X} \tfrac{\bar X}{\cancel{s_X}} \\
&= \alpha + \cancel{\beta \bar X} + \beta X - \cancel{\beta \bar X}
\end{align}$$
You can see the same with a data example:
> set.seed(42)
> n <- 200
> x1 <- rnorm(n, 1, 20)
> x2 <- rnorm(n, 2, 40)
> y <- 5 + 3*x1 + 2*x2 + rnorm(n, 0, 30)
> lm(y ~ x1 + x2)

Call:
lm(formula = y ~ x1 + x2)

Coefficients:
(Intercept)           x1           x2  
      3.180        2.922        2.050  

> lm(y ~ I((x1-1)/20) + I((x2-2)/40))

Call:
lm(formula = y ~ I((x1 - 1)/20) + I((x2 - 2)/40))

Coefficients:
   (Intercept)  I((x1 - 1)/20)  I((x2 - 2)/40)  
         10.20           58.44           82.02  

> 3.17952 + 58.4442*(1/20) + 82.015*(2/40)
[1] 10.20248
> 2.92221 * 20
[1] 58.4442
> 2.05038 * 40
[1] 82.0152

So as you can see:

*

*After the scaling, the parameters can be higher than three, or lower than negative three.

*The values of the parameters would depend on the values of the parameters before scaling and on the standard deviations of the raw data.

There is really no reason for assuming such bounds.
