> ## Documentation Index
> Fetch the complete documentation index at: https://kstats.oremif.org/llms.txt
> Use this file to discover all available pages before exploring further.

# Descriptive Statistics

> Summary statistics, central tendency, dispersion, quantiles, shape measures, frequency analysis, and streaming calculations in kstats-core.

{/*---IMPORT org.oremif.kstats.sampling.samples.DocsSamples--*/}

`kstats-core` provides descriptive statistics as extension functions on `DoubleArray` and `Iterable<Double>`. It also includes `Int` and `Long` overloads for the most common operations.

## Summary Snapshot

`describe()` is the fastest way to get a complete picture of a sample. It returns a `DescriptiveStatistics` object combining count, center, spread, quartiles, shape, and standard error.

```kotlin theme={"system"}
val data = doubleArrayOf(2.0, 4.0, 4.0, 4.0, 5.0, 5.0, 7.0, 9.0)
val stats = data.describe()

stats.count              // 8
stats.mean               // 5.0
stats.median             // 4.5
stats.standardDeviation  // 2.1380
stats.variance           // 4.5714
stats.q1                 // 4.0
stats.q3                 // 5.5
stats.interquartileRange // 1.5
stats.skewness           // 0.6563
stats.kurtosis           // -0.1640
stats.range              // 7.0
stats.standardError      // 0.7559
```

<Note>
  `variance()` and `standardDeviation()` use the **sample** formula (dividing by $n - 1$) by default. Pass `PopulationKind.POPULATION` to use the population formula (dividing by $n$).
</Note>

## Central Tendency

Different notions of "center" suit different data shapes. `mean()` is the arithmetic average, `median()` is the positional center, `mode()` returns the most frequent values, and the trimmed and weighted variants handle outliers and importance weights.

```kotlin theme={"system"}
val data = doubleArrayOf(2.0, 4.0, 4.0, 4.0, 5.0, 5.0, 7.0, 9.0)

data.mean()                    // 5.0
data.median()                  // 4.5
data.toList().mode()           // {4.0}

data.trimmedMean(0.1)          // 4.8333 — trims 10% from each tail

val positive = doubleArrayOf(1.0, 2.0, 4.0, 8.0)
positive.geometricMean()       // 2.8284
positive.harmonicMean()        // 2.1333

val values = doubleArrayOf(1.0, 2.0, 3.0)
val weights = doubleArrayOf(3.0, 1.0, 1.0)
values.weightedMean(weights)   // 1.6
```

* **Geometric mean** — use for growth rates and ratios. Requires all positive values.
* **Harmonic mean** — use for rates and speed averages. Requires all positive values.
* **Trimmed mean** — drops a proportion from each tail before averaging. Robust to outliers.
* **Weighted mean** — each observation contributes proportionally to its weight.

<Accordion title="Math details">
  $$
  \bar{x} = \frac{1}{n}\sum_{i=1}^{n} x_i
  $$

  $$
  \bar{x}_{\text{geo}} = \left(\prod_{i=1}^{n} x_i\right)^{1/n}
  $$

  $$
  \bar{x}_{\text{harm}} = \frac{n}{\sum_{i=1}^{n} \frac{1}{x_i}}
  $$

  $$
  \bar{x}_w = \frac{\sum_{i=1}^{n} w_i x_i}{\sum_{i=1}^{n} w_i}
  $$

  The trimmed mean drops $\lfloor p \cdot n \rfloor$ observations from each end of the sorted sample, then computes the arithmetic mean of the remaining values.
</Accordion>

## Dispersion

Spread measures describe how far values deviate from the center.

```kotlin theme={"system"}
val data = doubleArrayOf(2.0, 4.0, 4.0, 4.0, 5.0, 5.0, 7.0, 9.0)

data.variance()                                    // 4.5714 (sample)
data.variance(PopulationKind.POPULATION)            // 4.0
data.standardDeviation()                            // 2.1380
data.range()                                        // 7.0
data.interquartileRange()                           // 1.5
data.meanAbsoluteDeviation()                        // 1.5

data.trimmedVariance(0.1)                           // trimmed sample variance

data.semiVariance(5.0, SemiVarianceDirection.DOWNSIDE) // 1.7143
data.semiVariance(5.0, SemiVarianceDirection.UPSIDE)   // 2.8571
```

<Note>
  `semiVariance` measures variability on one side of a threshold. `DOWNSIDE` captures risk below the threshold; `UPSIDE` captures variability above it.
</Note>

<Accordion title="Math details">
  $$
  s^2 = \frac{1}{n-1}\sum_{i=1}^{n}(x_i - \bar{x})^2, \qquad
  s = \sqrt{s^2}, \qquad
  SE = \frac{s}{\sqrt{n}}
  $$

  $$
  \text{MAD} = \frac{1}{n}\sum_{i=1}^{n}|x_i - \bar{x}|
  $$

  $$
  \text{SemiVar}_{\text{down}}(t) = \frac{1}{n-1}\sum_{x_i < t}(x_i - t)^2
  $$
</Accordion>

## Quantiles and Position

Quantiles split the data at specified probability thresholds. `percentile()` takes a value on the 0–100 scale; `quantile()` takes a value on the 0–1 scale.

```kotlin theme={"system"}
val data = doubleArrayOf(1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0)

data.quantile(0.5)         // 5.5 (median)
data.quantile(0.25)        // 3.25 (Q1)
data.percentile(90.0)      // 9.1 (90th percentile)

val (q1, median, q3) = data.quartiles()
// q1 = 3.25, median = 5.5, q3 = 7.75
```

<Tip>
  `quantile()` supports a `QuantileInterpolation` parameter with options: `LINEAR` (default), `LOWER`, `HIGHER`, `NEAREST`, and `MIDPOINT`. These control how the value is interpolated when the quantile falls between two data points.
</Tip>

<Accordion title="Math details">
  With linear interpolation (the default), the quantile at probability $p$ for sorted data $x_{(1)}, \ldots, x_{(n)}$ is computed as:

  $$
  Q(p) = x_{(\lfloor h \rfloor)} + (h - \lfloor h \rfloor)(x_{(\lceil h \rceil)} - x_{(\lfloor h \rfloor)})
  $$

  where $h = p(n - 1) + 1$.
</Accordion>

## Shape and Moments

Shape measures describe asymmetry and tail weight beyond what the mean and variance capture.

```kotlin theme={"system"}
val data = doubleArrayOf(2.0, 4.0, 4.0, 4.0, 5.0, 5.0, 7.0, 9.0)

data.skewness()            // 0.6563 — positive: right tail is longer
data.kurtosis()            // -0.1640 — negative excess: lighter tails than normal
data.kurtosis(excess = false) // raw kurtosis (not centered at 0)

data.centralMoment(2)      // 4.0 (equals population variance)
data.centralMoment(3)      // 5.25
data.centralMoment(4)      // 44.5

data.kStatistic(1)          // 5.0 (equals mean)
data.kStatistic(2)          // 4.5714 (equals sample variance)
```

Skewness near zero suggests symmetry. Positive skewness indicates a longer right tail; negative indicates a longer left tail. Excess kurtosis near zero is similar to a normal distribution; positive values indicate heavier tails.

<Warning>
  `kurtosis()` returns **excess kurtosis** by default (`excess = true`), meaning the value is shifted so that a normal distribution has excess kurtosis of 0. Pass `excess = false` for the raw fourth standardized moment.
</Warning>

<Accordion title="Math details">
  $$
  g_1 = \frac{n}{(n-1)(n-2)} \sum_{i=1}^{n} \left(\frac{x_i - \bar{x}}{s}\right)^3
  $$

  $$
  g_2 = \frac{n(n+1)}{(n-1)(n-2)(n-3)} \sum_{i=1}^{n} \left(\frac{x_i - \bar{x}}{s}\right)^4 - \frac{3(n-1)^2}{(n-2)(n-3)}
  $$

  The $r$-th central moment is $\mu_r = \frac{1}{n}\sum_{i=1}^{n}(x_i - \bar{x})^r$. K-statistics ($k_r$) are unbiased estimators of the corresponding cumulants.
</Accordion>

## Frequency Analysis

The `Frequency` class counts occurrences, proportions, and cumulative frequencies for any `Comparable` type. `frequencyTable()` groups numeric data into equal-width bins.

```kotlin theme={"system"}
val freq = listOf("a", "a", "b", "b", "b", "c").toFrequency()
freq.totalCount            // 6
freq.count("b")            // 3
freq.proportion("b")       // 0.5
freq.cumulativeCount("b")  // 5 (a=2 + b=3)
freq.mode                  // {b}
```

{/*---FUN coreFrequencyTable--*/}

```kotlin theme={"system"}
val data = doubleArrayOf(1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0)
val bins = data.asIterable().frequencyTable(3)
// Each FrequencyBin has: range, count, relativeFrequency, cumulativeFrequency
bins[0].count              // number of values in the first bin
bins[0].relativeFrequency  // proportion of total
```

<Tip>
  `Frequency` works with any `Comparable` type — strings, enums, integers, or custom classes. For numeric binning, use `frequencyTable()` with either a bin count or a bin width.
</Tip>

## Streaming Statistics

`OnlineStatistics` computes mean, variance, skewness, and kurtosis incrementally without storing the full sample. Values can be added one at a time or in batches.

```kotlin theme={"system"}
val online = OnlineStatistics()
online.addAll(doubleArrayOf(2.0, 4.0, 4.0, 4.0, 5.0, 5.0, 7.0, 9.0))

online.count               // 8
online.mean                // 5.0
online.sum                 // 40.0
online.min                 // 2.0
online.max                 // 9.0
online.variance()          // 4.5714 (sample)
online.standardDeviation() // 2.1380
online.skewness()          // 0.6563

// Add more data later
online.add(3.0)
online.count               // 9
online.mean                // 4.7778
```

Use `OnlineStatistics` when data arrives incrementally (streaming, event-driven systems) or when buffering the full sample in memory is impractical.

<Accordion title="Math details">
  `OnlineStatistics` uses Welford's algorithm with the Terriberry extension for numerically stable single-pass computation of higher-order moments. The algorithm updates running sums of powers of deviations from the current mean, avoiding the catastrophic cancellation that affects naive two-pass formulas on large datasets.
</Accordion>

## Process Capability

Process capability indices quantify how well a process fits within its <Tooltip tip="The lower (LSL) and upper (USL) spec limits bound the range of acceptable values for the measured characteristic">specification limits</Tooltip>. `processCapability(lsl, usl)` returns four Statistical Process Control (SPC) indices in one call:

* **Cp, Cpk** — potential and actual capability using the sample standard deviation (short-term, within-subgroup spread).
* **Pp, Ppk** — the same formulas using the population standard deviation (long-term, overall spread).

The `-k` variants penalize a process that is off-center relative to the midpoint of the spec, so `Cpk` is never larger than `Cp`.

```kotlin theme={"system"}
// Ten parts measured against a spec window of [48, 52]
val measurements = doubleArrayOf(
    50.0, 50.5, 49.5, 50.2, 49.8, 50.1, 49.9, 50.3, 49.7, 50.0
)
val capability = measurements.processCapability(lsl = 48.0, usl = 52.0)

capability.cp   // 2.2646 — potential capability (spread vs tolerance)
capability.cpk  // 2.2646 — actual capability (penalizes off-centering)
capability.pp   // 2.3870 — overall (population σ) counterpart of Cp
capability.ppk  // 2.3870 — overall counterpart of Cpk
```

<Note>
  Use these indices only on a process that is already in statistical control (stable over time — see [Shewhart Control Charts](#shewhart-control-charts) below). For an unstable process the measured spread is not a fixed property of the process.
</Note>

Values `≥ 1.33` are usually considered capable; `≥ 1.67` highly capable. When `Cpk ≪ Cp`, re-center the process before trying to reduce variance.

<Accordion title="Math details">
  $$
  \mathrm{Cp} = \frac{\mathrm{USL} - \mathrm{LSL}}{6\sigma_s}, \qquad
  \mathrm{Cpk} = \min\!\left(\frac{\mathrm{USL} - \bar{x}}{3\sigma_s},\ \frac{\bar{x} - \mathrm{LSL}}{3\sigma_s}\right)
  $$

  Pp and Ppk use the population standard deviation $\sigma_p$ (divisor $n$) in place of the sample standard deviation $\sigma_s$ (divisor $n-1$). `processCapability` computes both in a single numerically stable Welford pass.
</Accordion>

## Shewhart Control Charts

Shewhart <Tooltip tip="Time-series plot of subgroup statistics with ±3σ control limits around a center line">control charts</Tooltip> plot subgroup statistics over time with three-sigma control limits. `xBarRChart()` monitors the process mean together with the range within each subgroup; `xBarSChart()` uses the sample standard deviation instead — more efficient for subgroup sizes above 10. Both require equal-sized subgroups of 2–25 observations.

```kotlin theme={"system"}
// Five subgroups of four parts; bracket width monitored per batch
val subgroups = listOf(
    doubleArrayOf(72.0, 84.0, 79.0, 49.0),
    doubleArrayOf(56.0, 87.0, 33.0, 42.0),
    doubleArrayOf(55.0, 73.0, 22.0, 60.0),
    doubleArrayOf(44.0, 80.0, 54.0, 74.0),
    doubleArrayOf(97.0, 26.0, 48.0, 58.0),
)
val chart = xBarRChart(subgroups)

chart.centerLine     // 59.65 — grand mean (x-double-bar)
chart.ucl            // 95.6626 — upper control limit for the mean
chart.lcl            // 23.6374 — lower control limit for the mean
chart.rChart.centerLine // 49.4 — average range (R-bar)
chart.rChart.ucl     // 112.7308 — upper limit for within-subgroup range
chart.rChart.lcl     // 0.0 — lower limit (D₃ = 0 for n ≤ 6)
```

Control limits use the standard SPC constants $A_2, A_3, D_3, D_4, B_3, B_4, c_4$ (Montgomery, *Introduction to Statistical Quality Control*, Appendix VI) available for subgroup sizes 2–25 through `spcConstants(n)`.

<Accordion title="Math details">
  For $k$ subgroups of size $n$ with subgroup means $\bar{x}_i$, ranges $R_i$, and standard deviations $s_i$:

  $$
  \text{x̄-R:}\quad \mathrm{UCL}/\mathrm{LCL} = \bar{\bar{x}} \pm A_2 \bar{R}, \quad R\text{-chart: } [D_3 \bar{R},\ D_4 \bar{R}]
  $$

  $$
  \text{x̄-S:}\quad \mathrm{UCL}/\mathrm{LCL} = \bar{\bar{x}} \pm A_3 \bar{s}, \quad S\text{-chart: } [B_3 \bar{s},\ B_4 \bar{s}]
  $$
</Accordion>

## CUSUM Chart

A Shewhart chart reacts slowly to drifts of less than 2σ because every point is judged in isolation. `cusum()` accumulates deviations from target over time, so a 0.5σ–1σ drift is detected within a few observations. The two-sided tabular form tracks an upper sum $C^+$ that catches upward shifts and a lower sum $C^-$ for downward shifts. An alarm fires on the first index where either sum exceeds the decision interval $H$.

```kotlin theme={"system"}
// Individual measurements from a process with target 10, drifting upward
val observations = doubleArrayOf(10.2, 10.4, 10.6, 10.9, 11.2, 11.5, 11.8, 12.0)
val result = cusum(observations, target = 10.0, k = 0.5, h = 3.0)

result.sPlus      // [0.0, 0.0, 0.1, 0.5, 1.2, 2.2, 3.5, 5.0]
result.sMinus     // all zero — no downward drift
result.alarmIndex // 6 — first index where C⁺ > H
```

<Tip>
  Tune `k` to half the shift size you want to detect, in units of σ — a common default is $K \approx 0.5\sigma$ targeting a 1σ drift. Set `h` to 4σ–5σ to match the in-control run length of a 3σ Shewhart chart while reacting much faster to small shifts.
</Tip>

<Accordion title="Math details">
  $$
  C^+_i = \max\bigl(0,\; C^+_{i-1} + (x_i - \mu_0 - K)\bigr), \qquad
  C^-_i = \max\bigl(0,\; C^-_{i-1} + (\mu_0 - K - x_i)\bigr)
  $$

  starting from $C^\pm_0 = 0$, with an alarm the first time $C^+_i > H$ or $C^-_i > H$.
</Accordion>

## EWMA Chart

EWMA (Roberts, 1959) is the other classic small-shift detector. Instead of an unbounded cumulative sum, `ewma()` maintains a weighted moving average that gives recent observations more influence while retaining memory of the past. Control limits widen with time until they reach a steady state, so the chart is most sensitive early — useful for catching an initial shift.

```kotlin theme={"system"}
// EWMA chart: target = 25, σ = 1, λ = 0.2, L = 3
val observations = doubleArrayOf(25.0, 24.5, 25.2, 26.1, 25.8, 27.0, 26.5, 28.0)
val result = ewma(
    observations,
    target = 25.0,
    sigma = 1.0,
    lambda = 0.2,
    controlLimitWidth = 3.0
)

result.smoothedValues[0] // 25.0 — Z₀ = λ·x + (1-λ)·target
result.smoothedValues[7] // 26.2549 — smoothed statistic at t = 7
result.ucl[0]            // 25.6 — narrow at first, widens with t
result.ucl[7]            // 25.9858 — approaching steady state
result.outOfControl      // [7] — Z₇ exceeds UCL₇
```

<Tip>
  $\lambda = 0.2$ with $L \approx 2.7$–$3.0$ is a common default. Smaller $\lambda$ emphasizes memory and detects smaller shifts; $\lambda = 1$ collapses EWMA into a Shewhart individuals chart.
</Tip>

<Accordion title="Math details">
  $$
  Z_t = \lambda x_t + (1 - \lambda) Z_{t-1}, \qquad Z_0 = \mu_0
  $$

  $$
  \mathrm{UCL}_t/\mathrm{LCL}_t = \mu_0 \pm L \sigma \sqrt{\frac{\lambda}{2 - \lambda}\,\bigl(1 - (1 - \lambda)^{2t}\bigr)}
  $$
</Accordion>

## Western Electric Rules

`westernElectricRules()` extends a Shewhart chart beyond the basic ±3σ check with four run-length heuristics that catch trends, clusters, and prolonged one-sided runs.

| Rule  | Pattern                                             | Detects                        |
| ----- | --------------------------------------------------- | ------------------------------ |
| **1** | 1 point beyond $\pm 3\sigma$                        | extreme single excursion       |
| **2** | 2 of last 3 points beyond $\pm 2\sigma$, same side  | strong shift                   |
| **3** | 4 of last 5 points beyond $\pm 1\sigma$, same side  | moderate shift                 |
| **4** | 8 consecutive points on the same side of the center | sustained shift, any magnitude |

Each rule's array contains the *trigger indices* — the observation whose arrival completes the offending pattern.

```kotlin theme={"system"}
// Process drifting upward in the last four observations
val observations = doubleArrayOf(
    0.1, 0.2, -0.3, 0.0, 1.4, 1.2, 2.4, 2.6, 3.5, 2.2
)
val violations = westernElectricRules(observations, center = 0.0, sigma = 1.0)

violations.rule1 // indices of points beyond ±3σ
violations.rule2 // indices where 2 of last 3 points are beyond ±2σ (same side)
violations.rule3 // indices where 4 of last 5 points are beyond ±1σ (same side)
violations.rule4 // indices where 8 consecutive points fall on the same side
```

<Tip>
  Combine a Shewhart chart (large shifts) with CUSUM or EWMA (small shifts) and Western Electric Rules (patterns) — the three views together catch the widest range of out-of-control conditions.
</Tip>

## Error Handling

<Warning>
  Empty arrays throw `InsufficientDataException`. Functions that require at least two observations (variance, standard deviation, skewness) throw `InsufficientDataException` for single-element input. Functions requiring positive data (`geometricMean`, `harmonicMean`) throw `InvalidParameterException` for non-positive values.
</Warning>

## API Reference

<Card title="Full API Reference" icon="book-open" href="https://oremif.github.io/kstats">
  Browse all public types, functions, parameter overloads, and return types in the Dokka-generated reference.
</Card>