Is the normality of residuals necessary to accept the null model in a multiple regression analysis? Best-subset regression analysis:
I want to test effects of differents ecological variables on my response variable. I am working with function glmulti() of glmulti R-package (method=gaussian and based on aicc values). When I run my scrip, the result output suggest mes the null model (y~1) as best model. 
I know that normality of residuals is important in regression analyses to accept the testing model but I'm not sure if you need normality too when your best model is the null model. My results suggest that none of my studied variable have effects on my Y variable, so...I think that it have sense that residuales of null model don't adjust with a normal distribution, what do you think?
This is part of my output with my 3 valid models (by differences <2 unids of AICc among themselves):
    model        aicc       weights
[1] Y ~ 1        33.79708   1.47E-01 #Best model = Null model
[2] Y ~ 1 + F2   34.84813   8.69E-02
[3] Y ~ 1 + F1   35.44111   6.46E-02

Model 1 (Y ~ 1) and normality test (below):
Coefficients:                   
            Estimate  Std. Error  t value   Pr(>|t|)    
(Intercept) 8.9914    0.1083      83.05     <2e-16  ***

---                 
Shapiro-Wilk normality test
data:  residuals((test1@objects[[1]]))
W = 0.78784, p-value = 0.0004302  --> No normality

Model 2 (Y ~ 1 + F2) and normality test (below):
Coefficients:           
            Estimate Std. Error t value Pr(>|t|)            
(Intercept)   8.9914     0.1067  84.283   <2e-16 ***
F2           -0.1381     0.1093  -1.263    0.222

---         
Shapiro-Wilk normality test
data:  residuals((test1@objects[[2]]))          
W = 0.93545, p-value = 0.1769  --> Normality    

Model 3 (Y ~ 1 + F2) and normality test (below):
Coefficients:               
            Estimate    Std. Error  t value   Pr(>|t|)
(Intercept) 8.9914      0.1082      83.102    <2e-16    ***
F1         -0.1121      0.1109      -1.011    0.325

---                 
Shapiro-Wilk normality test         
data:  residuals((test1@objects[[3]]))          
W = 0.86405, p-value = 0.007489 --> No normality

Note: My second model have normal residuals, but  F2 is not significant when test the model (only Intercept is significant). This is the same conclusion that the first model show, although this first has no normal residuals. Model 3 is not better than Model 2 because AICc is bigger and also show no-normal residuals.
By this reason I think that I could accept the null model (first model) as best model, what do you think?
Thanks for all.
 A: As a general rule, goodness-of-fit tests comparing regression models are robust to non-normality of the underlying error terms.  The reason for this is that most goodness-of-fit statistics are summation statistics that are subject to the central-limit-theorem.  It is not necessary for the underlying errors to be normally distributed for the goodness-of-fit statistic to converge in distribution to the distribution used for the test.
For a regression with an intercept and $k$ explanatory variables the AICc statistic can be written as:
$$\begin{equation} \begin{aligned}
\text{AICc} 
&= \frac{2nk}{n-k-1} - 2 \hat{\ell}_\mathbb{x,y} \\[6pt]
&= \frac{2nk}{n-k-1} - 2 \sum_{i=1}^n \ln p(\mathbf{x}_i, y_i | \hat{\boldsymbol{\beta}}, \hat{\sigma}) \\[6pt]
&= \frac{2nk}{n-k-1} + n \ln(2 \pi) + 2n \ln(\hat{\sigma}) + \frac{1}{\hat{\sigma}^2} \underbrace{\sum_{i=1}^n ( y_i - \mathbf{x}_i \cdot \hat{\boldsymbol{\beta}} )^2}_{\text{SSE}}. \\[6pt]
\end{aligned} \end{equation}$$
You can see from this expression that the AICc involves the residual-sum-of-squares (SSE).  For large $n$ the estimators of parameters converge to their true values, so they are not affected much by individual observations.  In this case, the SSE is a sum of (almost) IID random variables, and under some broad assumptions, this converges in distribution to the chi-squared distribution (which converges to normal).  Further information on the distribution of the AIC statistic can be found in Yanagihara and Ohmoto (2005), and this can easily be adjusted to the AICc.
Consequently, we can obtain a reasonable distributional approximation for the AICc statistic, even if the underlying error terms are not normally distributed.  We still require some conditions on the underlying error terms so that the central limit theorem applies; for example, we will generally require these to have finite variance, which rules out cases of heavy-tailed error distributions.
The above means that your model comparisons are probably quite robust to non-normality of the error terms, so long as $n$ is not too small, and so long as the broad conditions for applying the CLT are present.  Of course, regardless of which model you end up using, if the residuals show substantial departure from normality, then this means that the normality assumption in the model is false.  You might be able to improve your model using a GLM with a different error distribution, but even without this, many aspects of the model are robust.
