As far as I understand the sum(variances) should always be equal to the sum(explained_variance_), however when I run:
data = [[1, 4], [5, 1], [1, 4], [6, 8], [7, 1], [2,3], [3, 4], [1, 5], [3, 9]]
model = PCA(2)
model.fit(data)
variances = np.var(data, axis=0)
print(sum(variances), sum(model.explained_variance_))
I get this output:
11.283950617283951 12.694444444444443
I'm probably making some silly mistake but I don't see where is the problem.
sklearn.decomposition.PCA assumes that the model applied to a sample taken from a population, therefore when estimating the variances (actually the covariance matrice, but the variance sits in the diagonal) for unbiased estimation it uses $n - 1$ instead of $n$ as a divisor.
On the other hand, numpy.var with the default parameters calculates the population variance (with a divisor $n$), hence applied to a sample it's a biased estimation of the population variance.
To overcome this you have to use
variances = np.var(data, axis=0, ddof=1)
which is the same as
n= len(data)
variances = np.var(data, axis=0) * n / (n - 1)
In case if it's not a sample, but a full population (which is not a common use case), you have to amend the variances provided by the PCA model to be population variance by multiplying them with $\frac{n - 1}{n}$:
n= len(data)
sum(model.explained_variance_ * (n - 1) / n)
For an unknown reason, you have to multiply the sklearn explained_variance_ by $\frac{n - 1}{n}$ where $n$ is the number of observations to get the actual variance value. In your case, you have $n = 9$ observations in your data so $12.69 \cdot \frac{8}{9} = 11.28$
