# Expression for variation in the fourth moment of standard normal vs. sample size

Question:
What is the expression for the variation in the fourth (not-centralized) moment as a function of sample size for the standard normal distribution and is there a proper analytic/symbolic way to derive it?

Details:
When I use random number generator to create samples from a standard normal distribution, then I compute the mean, and I do this many times per sample-level, I can create a plot like the following:

My textbooks told me that the relationship between the variation and the sample size was:

$$Err \propto \frac{1}{\sqrt{n}}$$

where $Err$ is the standard deviation of the estimate and $n$ is the sample size.

When I plug it into this, I get a visual confirmation.

When I plug into a linear fit I get the following summary:

> summary(est)

Call:
lm(formula = s2 ~ I(sqrt(1/n)))

Residuals:
Min         1Q     Median         3Q        Max
-0.0097079 -0.0004087  0.0000041  0.0004380  0.0093119

Coefficients:
Estimate Std. Error  t value Pr(>|t|)
(Intercept)  -2.633e-04  5.312e-05   -4.957  8.4e-07 ***
I(sqrt(1/n))  1.005e+00  7.212e-04 1393.680  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.0009955 on 994 degrees of freedom
Multiple R-squared:  0.9995,    Adjusted R-squared:  0.9995
F-statistic: 1.942e+06 on 1 and 994 DF,  p-value: < 2.2e-16


The R-squared of 99.95% is an indicator of a pretty good fit.

What is the expression for the variation in the fourth (not-centralized) moment as a function of sample size for the standard normal distribution?

The plot of the variation in fourth moment versus sample size looks like this:

When I plot Error versus square root of inverse sample count, I do not get a straight line.

I could wire a stack of transformations on "n" into "glmulti" (link) but that can be a lot of work for little return.

Some "hacking" in algebra gives me a fair error expression, but I have no way of telling if it is local to my data, or more globally valid.

Here is my source code:

set.seed(5) #for reproducibility

n <- seq(from=5,to=1000,by = 1)

#number of iterations at each sample size
N <- 3000

#predeclare
store1 <- matrix(0, nrow=length(n),ncol=N)
store2 <- matrix(0, nrow=length(n),ncol=N)
mymean <- numeric()
mymom<- numeric()

nn <- matrix(0, nrow=length(n),ncol=N)

#for each sample size
for(i in 1:length(n)){

#repeat the measure "N" times
for (j in 1:N){

#store so we can separate it out later
nn[i,j] <- n[i]

#take sample
y <- rnorm(n = nn[i,j],mean = 0,sd = 1)

#compute moment
mymean <- mean(y)
mymom <- mean( (y-mean(y))^4)

#store
store1[i,j] <- mymean
store2[i,j] <- mymom
}
}

##compute variation

#predeclare
s1 <- numeric()
s2 <- numeric()

s3 <- numeric()
s4 <- numeric()

#loop
for (i in 1:length(n)){

#find elements which have a particular sample size
idx <- which(nn==n[i],arr.ind=T)

#compute statistics of first moment (aka mean)
s1[i] <- mean(store1[idx])
s2[i] <- sd(store1[idx])

#compute statistics of fourth moment (not centered
s3[i] <- mean(store2[idx])
s4[i] <- sd(store2[idx])

}

#make figures

plot(n, s2,xlab="Sample count",ylab="Error in mean")
lines(lowess(x=n,y=s2,f=0.01),col="Red",lwd=2)
grid()

plot(sqrt(1/n), s2,xlab="SQRT inverse Sample count",ylab="Error in mean")
lines(lowess(x=sqrt(1/n),y=s2,f=0.01),col="Red",lwd=2)
grid()

plot(n, s4,xlab="Sample count",ylab="Error in fourth moment")
lines(lowess(x=n,y=s4,f=0.005),col="Red",lwd=2)
grid()

plot(sqrt(1/n), s4,xlab="SQRT inverse Sample count",ylab="Error in fourth moment",
xlim=c(0,0.2),ylim=c(0,2))
lines(lowess(x=sqrt(1/n),y=s4,f=0.02),col="Red",lwd=2)
grid()

#fit to models
est <- lm(s2~I(sqrt(1/n)))
summary(est)

est2 <- lm(s2~I(sqrt(1/n))+I((1/n)^(1/4)) )
summary(est2)


## 1 Answer

What you are doing is basically correct, the small deviations you observe in the last plot is due to the large sample variation of empirical fourth moments.

Here is the algebra for the standard normal case. I am using formulas for moments from https://en.wikipedia.org/wiki/Normal_distribution#Moments Let $X_1, \dotsc, X_n$ be iid observations from a standard normal distribution. The fourth empirical moment (about zero) is given by $$\DeclareMathOperator{\E}{\mathbb{E}} \frac{1}{n} \sum_i X_i^4.$$Its expected value can be calculated: $$\E \left\{ \frac1{n} \sum_i X_i^4 \right\}= 3$$ and its variance: $$\DeclareMathOperator{\Var}{\mathbb{Var}} \Var\left\{ \frac1{n} \sum_i X_i^4 \right\} = (\frac1{n})^2\sum_i \Var X_i^4= \frac{96}{n}$$ and some details of the last calculation here: $$\Var X_i^4 = \E \left( X_i^4 -\E X_i^4 \right)^2 = \E \left( X_i^4-3\right)^2 = \E X_i^8 - 2\cdot 3 \E X_i^4 + 9 = \sigma^8\cdot 105 -6 \cdot 3 + 9= 96$$