IOAI Training: Learn AI by solving problems

The problem in brief

You receive one year of daily temperatures (365 numbers in a NumPy file) and must compute six descriptive statistics: the mean, the standard deviation, the variance, the kurtosis, the skewness, and the cumulative moving average (CMA), the average of all days up to each day. The answer is a CSV with one row per day: the day's temperature, its CMA, and the five global statistics repeated on every row. The statement also mentions a "pollution" field, but the file contains only temperatures.

The data are the official Kaggle file. Kaggle's grader is not public and we did not submit anything. Every value below was computed twice, with library calls and with the formulas written out by hand, and both agree to $10^{-9}$ .

How you are scored

The Evaluation page only says that the submission is judged on the correctness of the computed values, with $-1$ for any value you could not compute. It gives no formula and no tolerance. On the public leaderboard the best score is 1.1, reached by two entries, followed by 0.9. The scale is not explained.

First look at the data

The 365 values run from $-12.4$ °C to $58.5$ °C, with no missing values. A histogram looks very much like a bell curve (right panel below). It is in fact exactly one: the array equals 20 + 10 * np.random.randn(365) after np.random.seed(42), to the last digit. The "temperatures" are a normal sample with true mean 20 and true standard deviation 10. That tells us what the answers should be close to, which is a useful sanity check.

Temperatures, CMA and distribution

The key idea

Every statistic here is a moment of the data. Write $\bar x$ for the mean and $m_k$ for the $k$ -th central moment:

\bar x = \frac1n \sum_{i=1}^n x_i, \qquad m_k = \frac1n \sum_{i=1}^n (x_i - \bar x)^k .

Then

\text{variance} = m_2, \qquad \text{standard deviation} = \sqrt{m_2},

\text{skewness} = \frac{m_3}{m_2^{3/2}}, \qquad \text{excess kurtosis} = \frac{m_4}{m_2^{2}} - 3 .

Skewness and kurtosis are moments of the standardised data $z_i = (x_i - \bar x)/\sqrt{m_2}$ , so they have no unit and do not change if you convert to Fahrenheit. The $-3$ makes the kurtosis of a normal distribution equal to 0.

The only real difficulty is that several conventions exist, and the statement does not pick one:

quantity	divide by $n$ (library default)	other common convention
variance	89.6338 (`np.var`)	89.8801, divide by $n-1$ (`ddof=1`, `pandas`)
standard deviation	9.4675 (`np.std`)	9.4805, divide by $n-1$
kurtosis	0.5913, excess (`scipy.stats.kurtosis`)	3.5913 Pearson (`fisher=False`), 0.6161 bias-corrected
skewness	0.1587 (`scipy.stats.skew`)	0.1593 bias-corrected

We use the defaults of NumPy and SciPy (left column), which treat the 365 days as the whole population. Beware of mixing libraries: pandas.Series.std() divides by $n-1$ , while np.std divides by $n$ , and pandas.Series.kurt() is the bias-corrected version.

The CMA has a neat closed form. The CMA of day $k$ is $\frac{1}{k}\sum_{i \le k} x_i$ , which is just the running sum divided by the day number:

cma = np.cumsum(t) / np.arange(1, n + 1)

That is one pass over the data. Recomputing the average from scratch for each day would cost $O(n^2)$ operations, which is the trap the statement's remark on efficiency points to. The equivalent incremental update, useful when data arrive one day at a time, is

\text{CMA}_k = \text{CMA}_{k-1} + \frac{x_k - \text{CMA}_{k-1}}{k}.

The full solution

Load the array with np.load(..., allow_pickle=False).
Compute the mean, np.std, np.var, scipy.stats.kurtosis and scipy.stats.skew, and the CMA with a cumulative sum.
Recompute everything from the moment formulas and with a day-by-day loop for the CMA, and assert that the two agree.
Write the 365-row CSV in the requested column order.

Results

statistic	value
mean	20.0995 °C
standard deviation (divide by $n$ )	9.4675 °C
variance (divide by $n$ )	89.6338 °C²
excess kurtosis	0.5913
skewness	0.1587
CMA on day 1 / day 2 / day 365	24.967 / 21.792 / 20.0995

The last CMA equals the mean, as it must.

Do these numbers make sense? The mean and standard deviation are close to the true 20 and 10: the standard error of the mean is $10/\sqrt{365} \approx 0.52$ , and we are 0.1 away. For a normal sample of size $n$ , skewness and excess kurtosis have standard errors of about $\sqrt{6/n} = 0.128$ and $\sqrt{24/n} = 0.256$ . Our skewness (0.16) is unremarkable. The kurtosis (0.59) is 2.3 standard errors above 0, and one single day explains it: the 58.5 °C day lies 4.06 standard deviations above the mean, and without it the kurtosis falls to 0.12 and the skewness to 0.01. The fourth power makes kurtosis extremely sensitive to the most extreme points. The left panel shows the CMA wandering early and settling as days accumulate, inside a band that shrinks like $\sigma/\sqrt{k}$ .

What to remember

Mean, variance, skewness and kurtosis are the first four moments, the last two computed on standardised data. Writing them by hand once is the best way to remember what each measures.
State your conventions. $n$ versus $n-1$ , excess versus Pearson kurtosis, biased versus corrected estimators: libraries disagree (NumPy divides by $n$ , pandas by $n-1$ ), and a grader expects one of them.
Cumulative quantities are cumulative sums. np.cumsum turns an $O(n^2)$ loop into $O(n)$ , and the incremental update $\text{CMA}_k = \text{CMA}_{k-1} + (x_k - \text{CMA}_{k-1})/k$ works on streaming data.
Higher moments are fragile. One extreme value out of 365 moved the kurtosis from 0.12 to 0.59. Look at the extremes before interpreting a kurtosis.
Sanity-check against what the data should be. Knowing the expected value and the standard error of each statistic tells you at once whether a result is plausible or a bug.

"""BdAIO 2024, "Earth's Pulse": worked solution.

Computes, for the 365 daily temperatures, the mean, standard deviation, variance, kurtosis, skewness
and the cumulative moving average, writes them in the submission format, and cross-checks every
statistic two ways (library call and the formula written out by hand).

Conventions (the statement does not fix them, we use the library defaults):
  standard deviation and variance: population versions (divide by n), as numpy does;
  kurtosis: EXCESS kurtosis (normal distribution = 0), biased estimator, as scipy.stats.kurtosis does;
  skewness: biased estimator, as scipy.stats.skew does.
The other conventions are printed too.

Usage: python solution.py [path/to/temperatures.npy]
"""
import sys
import numpy as np
import pandas as pd
from scipy import stats

path = sys.argv[1] if len(sys.argv) > 1 else "work/temperatures.npy"
t = np.load(path, allow_pickle=False)               # a plain float64 array: no pickle involved
n = len(t)

# ---------------------------------------------------------------- library versions
mean, std, var = t.mean(), t.std(), t.var()
kurt, skew = stats.kurtosis(t), stats.skew(t)
cma = np.cumsum(t) / np.arange(1, n + 1)           # CMA_k = (x_1 + ... + x_k) / k, in O(n)

# ---------------------------------------------------------------- the same, by hand
m = sum(t) / n
d = t - m
m2, m3, m4 = (d ** 2).mean(), (d ** 3).mean(), (d ** 4).mean()   # central moments
hand = {"mean": m, "var": m2, "std": m2 ** 0.5, "skew": m3 / m2 ** 1.5, "kurt": m4 / m2 ** 2 - 3}
lib = {"mean": mean, "var": var, "std": std, "skew": skew, "kurt": kurt}
for k in hand:
    assert abs(hand[k] - lib[k]) < 1e-9, k
running, cma_loop = 0.0, []
for k, x in enumerate(t, 1):                         # the definition, one day at a time
    running += x
    cma_loop.append(running / k)
assert np.allclose(cma_loop, cma)

print(f"mean {mean:.4f}  std {std:.4f}  var {var:.4f}  kurtosis {kurt:.4f}  skewness {skew:.4f}")
print(f"other conventions: std (n-1) {t.std(ddof=1):.4f}, var (n-1) {t.var(ddof=1):.4f}, "
      f"Pearson kurtosis {stats.kurtosis(t, fisher=False):.4f}, "
      f"unbiased kurtosis {stats.kurtosis(t, bias=False):.4f}, unbiased skewness {stats.skew(t, bias=False):.4f}")
print(f"standard errors under normality: skewness ~ {np.sqrt(6 / n):.3f}, kurtosis ~ {np.sqrt(24 / n):.3f}")

# the data are exactly 20 + 10 * N(0, 1) drawn with numpy's legacy seed 42
np.random.seed(42)
print("equals 20 + 10 * np.random.randn(365) with seed 42:", np.allclose(t, 20 + 10 * np.random.randn(n)))

sub = pd.DataFrame({"Day": np.arange(1, n + 1), "Temperature": t, "Cumulative_Moving_Average": cma,
                    "Mean_Temperature": mean, "Standard_Deviation": std, "Variance": var,
                    "Kurtosis": kurt, "Skewness": skew})
sub.to_csv("submission.csv", index=False)
print("wrote submission.csv (not submitted)")