6
$\begingroup$

I am fairly new to hypothesis testing. I read that the F-test for equality of variance doesn't work well for distributions that aren't normal. Does this test work for truncated normal distributions for cases when the truncated distribution is symmetric?

$\endgroup$
3
  • 2
    $\begingroup$ The F-test relies on the fact that the standard deviation of a normally distributed sample follows a (scaled) chi-square distribution. This fact will not hold for truncated normals, and so the F-test assumptions are not satisfied. Perhaps you could look into a flexible testing procedure, like the Bootstrap. $\endgroup$
    – knrumsey
    Apr 20 at 18:52
  • $\begingroup$ Thank you, can you please suggest some good books which present hypothesis testing? Mathematical rigour will be much appreciated. $\endgroup$
    – junfan02
    Apr 20 at 20:02
  • 3
    $\begingroup$ Casella and Berger (Statistical Inference) is a widely used text which you might find helpful. There are more rigorous books but it's quite good. $\endgroup$
    – Glen_b
    Apr 21 at 2:06

2 Answers 2

7
$\begingroup$

The F-test relies on the fact that the variance of a normally distributed sample follows a (scaled) chi-squared distribution. This fact does not hold for truncated normals, and thus the F-test assumptions are not satisfied.

Edit (4/25/2024): Added a simulation study to compare methods, per the discussion in the comments (@whuber, @StephenKolassa, @SextusEmpiricus).


Testing for Equal Variance with the Bootstrap

The bootstrap algrorithm is a powerful and flexible approach for hypothesis testing in general. There are lots of questions and answers on this site and others about how it works, so I'll skip the details here. I'll show you how it works in R. Lets start by generating some data.

set.seed(192012)
x <- rnorm(500, 0, 1)       # Mean 0, sd 1
x <- x[which(abs(x) < 1)]   # Truncated on [-1, 1] interval

y <- rnorm(500, 0, 0.6)     # Mean 0, sd 0.6
y <- y[which(abs(y) < 1.5)] # Truncated on [-1.5, 1.5] interval

The variance of the $y$ samples is greater than that of the $x$ samples, but we question whether this difference is "statistically significant".

> var(x)
[1] 0.2887654
> var(y)
[1] 0.3382733

We can approximate the sampling distribution of the sample variances via the bootstrap:

B <- 10000
boot_x <- boot_y <- rep(NA, B)
for(i in 1:B){
  ind_x <- sample(length(x), length(x), replace=TRUE)
  boot_x[i] <- var(x[ind_x])
  
  ind_y <- sample(length(y), length(y), replace=TRUE)
  boot_y[i] <- var(y[ind_y])
}

$95\%$ confidence intervals can be obtained for each sample with the command:

> quantile(boot_x, c(0.025, 0.975))
     2.5%     97.5% 
0.2593101 0.3171042 
> quantile(boot_y, c(0.025, 0.975))
     2.5%     97.5% 
0.3019288 0.3750763 

The fact that these intervals overlap considerably suggests that we don't quite have enough evidence to claim a difference between the variances of the two samples (at the $5\%$ significance level). Plotting the bootstrap distributions is always a good idea:

dx <- density(boot_x)
dy <- density(boot_y)
plot(NULL, xlim=range(c(dx$x, dy$x)), ylim=range(c(dx$y, dy$y)),
     xlab="", ylab="Density", main="Bootstrap Distributions")
polygon(dx, border='white', col=adjustcolor("orange", alpha.f=0.5))
polygon(dy, border='white', col=adjustcolor("dodgerblue", alpha.f=0.5))
legend("topright", c("Bootstrap Distribution of S_x", "Bootstrap Distribution of S_y"),
       fill=c("orange", "dodgerblue"))

enter image description here


Textbook Recommendations

For a some books which present hypothesis testing, consider

  • Statistical Rethinking by Richard McElreath. A modern stat inference book which I really like. Has a lot of code, and discusses both Bayesian and Frequentist frameworks. Not particularly mathematical.
  • Most classical inference textbooks will cover this with some rigor. See Casella & Berger or Linear Statistical Inference by C.R. Rao.
  • For Bootstrapping, people seem to like Davison and Hinkley (1997), "Bootstrap Methods and Their Application."

A Simulation Study

There has been some relevant discussion in the comments (e.g., @whuber, @stephenkolassa). To partially address this discussion, we can conduct a simulation study to compare some different methods:

  1. The bootstrap test
  2. F-test
  3. The modified F-test proposed by @SextusEmpiricus. In the following, let $S_j^2$ be the sample variance and $N_j$ be the sample size for sample $j=1,2$.

To compute a formal p-value for the Bootstrap procedure, we will keep track of two counts: $C_1$ is the number of times that $S^2_1 < S_2^2$ and $C_2$ is the number of times that $S^2_2 < S_1^2$. In a one-sided test, only one of these counts will matter. For a two-sided test, we can define the p-value as $$p_\text{boot} = \frac{2}{B}\min\{C_1, C_2\},$$ where $B$ is the number of bootstrap samples.

For the F-test, we compute the test statistic $$F = \frac{S_1^2}{S_2^2}$$ and we compute the p-value via $$p_F = 2\min\{P_F(F|N_1-1, N_2-1), 1 - P_F(F|N_1-1, N_2-1),$$ where $P_F$ is the distribution function of an F-distribution with $N_1-1$ and $N_2-1$ degrees of freedom.

In @SextusEmpiricus' answer, a modification is proposed which can better capture the distribution of the test statistic under the null hypothesis. In the case of unequal sample sizes, the corrected p-value can be computed via the formula $$p^*_F = 2\min\{P_F(F|\nu_1, \nu_2), 1 - P_F(F|\nu_1, \nu_2)\},$$ where $\nu_j = 3(N_j-1)/\hat\kappa_j$ is the corrected degrees of freedom and $\hat\kappa_j$ is an estimate of the kurtosis.

For anybody that wants to play around with different cases, the R code I use to conduct the simulation study can be found here. I will focus on 3 cases for simplicity:

  1. Mild truncation, slightly different (true) variances
  2. Mild truncation, roughly equal (true) variances
  3. Major, asymmetric, truncation. Different variances.

Case 1 and 2

For both of these cases, we set \begin{align*} &n_1=100, \ \mu_1 = 0, \ \sigma_1 = 1 \\ &n_2=150, \ \mu_2 = 0, \ \sigma_2 = 0.6, \ a_2 = -1.5, \ b_2 = 1.5 \end{align*} The truncation interval for the first sample is set to $$\begin{cases} (-1, 1), & \text{Case 1} \\ (-1.07265, 1.07625), & \text{Case 2}. \end{cases}$$ We choose the values so that $$\begin{cases} \text{Var}(X_1) \approx 0.2911 < 0.3281 \approx \text{Var}(X_2), & \text{Case 1} \\ \text{Var}(X_1) \approx 0.3281 \approx \text{Var}(X_2), & \text{Case 2}. \end{cases}$$

We simulate $100$ different datasets using the settings described above, and we compute the p-value for each of the three methods. The p-values are plotted below (conveniently ordered for visualization). Recall that under the null hypothesis, we expect the p-value to be uniformly distributed; this is represented below by the grey "reference line".

In case 1, we see that the Bootstrap seems to be the most powerful of the three methods, followed by the modified F-test. sim study case 1

In case 2, we expect the p-values to roughly align with the reference curve. In this setting, all methods seem to perform similarly, but the bootstrap method could be said to be the most "calibrated". sim study case 2

Case 3

We conclude with a case where one of the truncated normal distributions is highly non-normal. We use the same settings as above, except for the truncation interval for the first sample becomes $(0.5, 2.5).$ Again, we see that the three methods perform similarly. sim study case 3

$\endgroup$
7
  • 6
    $\begingroup$ "...assumptions not satisfied" isn't a relevant criterion, because it is (extremely) rare for such specific mathematical assumptions to be satisfied in any application. What matters is how and to what degree the assumptions are violated. With a truncated Normal, it's likely the violation is mild and that, if anything, an F test would work better in this case than in most cases, depending on how much truncation occurs and how much data there are. $\endgroup$
    – whuber
    Apr 21 at 3:31
  • 2
    $\begingroup$ It depends on the degree of truncation. $\endgroup$
    – whuber
    Apr 21 at 5:03
  • 2
    $\begingroup$ +1, it's always good to have alternatives that rely on fewer assumptions. That said, I would be interested in how the F test actually performed in this case, where the truncation is still not all that strong. Maybe add that? $\endgroup$ Apr 21 at 5:48
  • 2
    $\begingroup$ @StephanKolassa A good idea; took me a while to get to it, but I added a fairly complete comparison. $\endgroup$
    – knrumsey
    Apr 25 at 19:17
  • 1
    $\begingroup$ @Glen_b Error has been fixed, thanks! $\endgroup$
    – knrumsey
    Apr 25 at 19:28
3
$\begingroup$

I read that the F-test for equality of variance doesn't work well for distributions that aren't normal.

The F-test gives an exact result when the distributions are normal.

When this is not the case then the F-test will still be an approximation (a good one if the sample sizes are large).

The F-test relates to a ratio of two sums

$$F = \frac{ \frac{1}{n-1}\sum_{i=1}^n (X_i-\bar{X})^2}{ \frac{1}{m-1} \sum_{i=1}^m (Y_i-\bar{Y})^2}$$

those sums will both approach a normal distribution for large sample numbers, and the F-distribution with large degrees of freedom is more like a ratio of normal distributions.

The trick is to figure out the mean and variance of those normal distributions in the numerator and denominator. We can relate this to the second and fourth central moments $m_2$ and $m_4$ of $X$ and $Y$

$$\begin{array}{} \text{Mean} &=& m_2 \\ \text{Variance} &=& (m_4 - m_2^2)/\nu \\ &=& \frac{\kappa -1}{\nu} m_2^2 \end{array}$$

Where $\kappa = m_4/m_2^2$ is the kurtosis.

In the case of a normal distribution these will be $m_2 = \sigma^2$ and $m_4 - m_2^2 = 2\sigma^4$

As a demonstration see for example $F(40,40)$ compared with a related normal ratio distribution:

comparing f distribution with normal ratio distribution

So more important than the distribution being normal is that kurtosis is approximately equal to three. If this is not the case then one might adjust for this and adjust the degrees of freedom $$\nu' = \frac{3-1}{\kappa-1}\nu$$

Below this is demonstrated for $X \sim U(0,3\sqrt{2})$ and $Y \sim t_6$ which are two distributions with the same variance but a large difference in kurtosis. If the samples are large (in the image we uses $n_1 = n_2 = 501$), then the adjusted F-test more closely resembles the simulated tests.

testing F-test with adjustment for degrees of freedom

In the image the result is not perfect, but it is also an extreme example. If we model two times a uniform distribution then the adjusted f-distribution is much better and already at lower sample sizes.

So, for a truncated normal distribution it seems that one can use an f-test, but one may need to adjust the degrees of freedom. I believe that this can be easier than computing a bootstrap distribution.


R-code for the image 1

#
# function to compute Hinkley's density formula
#
fw <- function(w,mu1,mu2,sig1,sig2,rho) {
  #several parameters
  aw <- sqrt(w^2/sig1^2 - 2*rho*w/(sig1*sig2) + 1/sig2^2)
  bw <- w*mu1/sig1^2 - rho*(mu1+mu2*w)/(sig1*sig2)+ mu2/sig2^2
  c <- mu1^2/sig1^2 - 2 * rho * mu1 * mu2 / (sig1*sig2) + mu2^2/sig2^2
  dw <- exp((bw^2 - c*aw^2)/(2*(1-rho^2)*aw^2))
  
  # output from Hinkley's density formula
  out <- (bw*dw / ( sqrt(2*pi) * sig1 * sig2 * aw^3)) * (pnorm(bw/aw/sqrt(1-rho^2),0,1) - pnorm(-bw/aw/sqrt(1-rho^2),0,1)) + 
    sqrt(1-rho^2)/(pi*sig1*sig2*aw^2) * exp(-c/(2*(1-rho^2)))
  
  out
}
fw <- Vectorize(fw)



x = seq(0,3,0.01)

nu_1 = 40
nu_2 = 40

mu_1 = 1
mu_2 = 1
sig1 = sqrt(2/nu_1)
sig2 = sqrt(2/nu_2)

plot(x,df(x,nu_1,nu_2), type  = "l", 
    main = "comparing F distribution \n with normal ratio distribution",
    ylab = "density", xlab = "F statistic", lwd = 1.5)
lines(x,fw(x,mu_1,mu_2,sig1,sig2,0), col = 3, lwd = 1.5, lty = 2)

legend(1.35,1.27, 
       c("F(40,40) distribution","normal ratio distribution"),
       lty = c(1,2),
       col = c(1,3),
       cex = 0.7, lwd = 1.5)

R-code for image 2

set.seed(1)
n1 = 501
n2 = 501

fstat = replicate(10^5, expr = {
  x = runif(n1,0,sqrt(18))   ### variance = 1/12 b^2 = 1.5 ### kurtosis = 9/5 
  y = rt(n2,6)               ### variance = 1.5            ### kurtosis = 6

  (sum((x-mean(x))^2)/(n1-1))/
  (sum((y-mean(y))^2)/(n2-1))
})

hist(fstat, freq = 0, breaks = seq(0,max(fstat+0.05),0.02), xlim = c(0.5,2),
     main = "histogram of f-statistic compared with \n different f distributions",
     cex.main = 1)
x = seq(0,3,0.001)
lines(x,df(x, (n1-1)*(3-1)/(9/5-1) ,  (n2-1)*(3-2)/(6-1) ))
lines(x,df(x, (n1-1) ,  (n2-1) ), col = 3)


legend(1.2,4, 
       c("standard distribution",
         "adjusted distribution"),
       lty = c(1,2),
       col = c(1,3),
       cex = 0.7, lwd = 1.5)
$\endgroup$
6
  • 1
    $\begingroup$ Disclaimer, I am just making this idea about using an adjustment $\nu' = \frac{2}{Kurt}\nu$ up from scratch. It is not rigourously studied by me. I can imagine however that this correction factor has been described before. $\endgroup$ Apr 23 at 14:32
  • $\begingroup$ +1, a very cool idea. I am currently working on a comparison which includes this method. Clarification: You do mean that the kurtosis should be approximately $3$, right (for a normal distribution)? If not, then I don't know where the number $2$ comes from. $\endgroup$
    – knrumsey
    Apr 25 at 18:18
  • $\begingroup$ @knrumsey what an embarrassing mistake. Yes it is supposed to relate to the kurtosis of the normal distribution which is 3 instead of 2. $\endgroup$ Apr 25 at 18:59
  • 1
    $\begingroup$ It is a bit annoying however that now, when I use two uniform distributions, the adjusted degrees of freedom cases do not match so perfect as before. I have to rework this with some pen and paper, and coffee. (Maybe I am missing a square root or something like that) $\endgroup$ Apr 25 at 19:13
  • 1
    $\begingroup$ @knrumsey the factor 2 was right after all. I must have accidentally stumbled upon it because I did not do it on purpose, but the simulations worked well with it. I have now edited the question with a better explanation of the relationship with the kurtosis. While improving the description of the derivation I noticed that there should be a $-1$ term added to the equation. $\endgroup$ Apr 25 at 19:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.