Source: Tamhane/Dunlop Ch. 2, Sections 4.1, 4.2, 4.3.4 | Sheather (2004)
Tags: probability density function estimator, PMF estimator, histogram, relative frequency histogram, density estimator, kernel density estimator, KDE, bandwidth, Gaussian kernel, bin width, Scott's rule, Freedman-Diaconis rule, Silverman rule, MISE, AMISE
For discrete distributions, the PMF is estimated by dividing observed frequencies by the sample size. For continuous distributions, there are two main approaches: the relative frequency histogram (which divides the data range into bins and counts) and the kernel density estimator (which places a smooth bump at each data point and sums them). The kernel density estimator produces a smooth curve and is generally superior, but its quality depends heavily on the choice of bandwidth.
Probability mass function (PMF)
For a discrete random variable, f(y_j) = Pr[Y = y_j]. Estimated by the proportion of observations equal to y_j.
Probability density function (PDF)
For a continuous random variable, the derivative of the CDF: f(y) = dF(y)/dy. Estimated by measuring the local concentration of data near y.
Local density
The relative concentration of data in a small interval centred at y. Formally: f-hat(y) = (fraction of data in the interval) / (length of interval).
Relative frequency histogram (density histogram)
A histogram where the vertical axis shows (relative frequency) / (class width), so the total area under the bars equals 1. This makes it a proper density estimator.
Frequency histogram
A histogram where the vertical axis shows raw counts. Not a proper density estimator because the area under the bars does not equal 1.
Bin width (h) / class width / mesh size
The width of each interval in a histogram. A critical choice: too wide and you lose detail, too narrow and you get noise.
Kernel density estimator (KDE)
A smooth, continuous density estimate: f-hat(y) = (1/(nh)) × Σ K((y - Yi)/h), where K is a kernel function and h is the bandwidth.
Kernel function K
A non-negative function that integrates to 1, is symmetric about zero, and is typically unimodal. Common choices: Gaussian, Epanechnikov, triangular, rectangular, cosine.
Bandwidth (h) / smoothing parameter
Controls the width of each kernel bump. Small h produces a wiggly, undersmoothed estimate; large h produces a flat, oversmoothed estimate.
MISE (mean integrated squared error)
The expected integrated squared difference between the estimated and true density. The standard criterion for evaluating density estimators.
AMISE (asymptotic MISE)
The large-sample approximation to MISE: AMISE = (1/(nh)) × R(K) + (h^4/4) × μ_2(K)^2 × R(f'').
For a discrete random variable Y with possible values y_1, y_2, ..., y_k:
Count the observed frequency: f-hat_j = number of Yi equal to y_j
The estimated PMF at y_j is: f-hat(y_j) = f-hat_j / n
At all other values: f-hat(y) = 0
Two standard plots for discrete data:
Bar graph: bars centred at each y_j with height proportional to f-hat(y_j)
Line plot (needle plot): vertical lines from the x-axis up to f-hat(y_j) at each y_j
Example from the handout: 1000 containers inspected for defects (0 to 10 per container). The estimated PMF at 0 defects is 402/1000 = 0.402, at 1 defect is 321/1000 = 0.321, and so on.
The PDF is defined as f(y) = lim (Δ→0) [F(y + Δ/2) - F(y - Δ/2)] / Δ.
For a small interval of width h centred at y:
h × f(y) ≈ probability that Y falls in (y - h/2, y + h/2)
The sample-based estimator replaces the probability with the sample proportion:
f-hat(y) = (number of Yi in (y - h/2, y + h/2)) / (n × h)
Set up a mesh of k equal-width bins covering the data range:
Bin endpoints: t_0 = Y_(1), t_k = Y_(n)
Bin width: h = (Y_(n) - Y_(1)) / k = range / k
n_j = number of observations in the j-th bin [t_(j-1), t_j)
R_j = n_j / n (relative frequency)
f-hat_j = R_j / h (density estimate in the j-th bin)
The density estimator is:
f-hat(y) = f-hat_j if t_(j-1) ≤ y < t_j
f-hat(y) = 0 otherwise
Plots of n_j, R_j, and f-hat_j vs y all have the same shape when bins are equal width, but only f-hat_j is a proper density estimator (non-negative, total area = 1).
Method | Formula | Notes |
|---|---|---|
Textbook default | k between 5 and 15 | Simple but arbitrary |
R default | h = range / (1 + log_2(n)) | Sturges' rule |
Scott's rule | h = 3.5 × s × n^(-1/3) | s = sample standard deviation; use |
Freedman-Diaconis rule | h = 2 × IQR × n^(-1/3) | IQR = interquartile range; use |
When bins have different widths h_j = t_j - t_(j-1):
f-hat_j = n_j / (n × h_j), where h_j is the width of the j-th bin
With unequal widths, plots of raw counts n_j or relative frequencies R_j are misleading because wider bins accumulate more counts simply by being wider. Only the density estimate f-hat_j = R_j / h_j gives a fair picture.
It measures local density only at the midpoint of each bin, then assigns that value to the entire bin, producing a piecewise constant function
It uses only the data within each bin to estimate f(y) at points in that bin, ignoring all other data
The kernel density estimator addresses both of these issues
The kernel density estimate of f(y) is:
f-hat(y) = (1/(nh)) × Σ K((y - Yi)/h), summing from i = 1 to n
Each observation Yi contributes a scaled copy of the kernel K, centred at Yi and spread by bandwidth h. The estimate at any point y is the sum of all n contributions.
Key properties of f-hat:
f-hat(y) ≥ 0 for all y (because K is non-negative)
The integral of f-hat over the real line equals 1 (because K integrates to 1)
So f-hat is a valid PDF
Place a small bump (kernel) centred at each data point
The bump is wider or narrower depending on the bandwidth h
Sum all n bumps and divide by n
Data points far from y contribute very little to f-hat(y); nearby points contribute a lot
Kernel | Formula | Support |
|---|---|---|
Gaussian | (1/√(2π)) × exp(-u²/2) | All real numbers |
Epanechnikov | (3/4)(1 - u²) for |u| ≤ 1 | [-1, 1] |
Triangular | (1 - |u|) for |u| ≤ 1 | [-1, 1] |
Rectangular (uniform) | 1/2 for |u| ≤ 1 | [-1, 1] |
Quartic (biweight) | (15/16)(1 - u²)² for |u| ≤ 1 | [-1, 1] |
Cosine | (π/4) × cos(πu/2) for |u| ≤ 1 | [-1, 1] |
The Gaussian kernel is the most popular default. The choice of kernel matters far less than the choice of bandwidth.
Bandwidth controls the bias-variance tradeoff:
Small h (undersmoothing): low bias, high variance, wiggly estimate with spurious peaks
Large h (oversmoothing): high bias, low variance, overly flat estimate that misses real features
The formal expressions for bias and variance of f-hat(y) are:
Bias{f-hat(y)} ≈ (h²/2) × μ_2(K) × f''(y)
Var{f-hat(y)} ≈ (1/(nh)) × R(K) × f(y)
where R(K) = ∫K²(y)dy and μ_2(K) = ∫u² K(u) du.
The overall measure of estimation quality is:
MISE(f-hat) = E[∫(f-hat(y) - f(y))² dy] = ∫Bias² dy + ∫Variance dy
The asymptotic version is:
AMISE(f-hat) = (1/(nh)) × R(K) + (h⁴/4) × μ_2(K)² × R(f'')
where R(f'') = ∫[f''(y)]² dy measures the roughness (curvature) of the true density. More curvature means harder estimation and a smaller optimal bandwidth.
h_AMISE = [R(K) / (μ_2(K)² × R(f''))]^(1/5) × n^(-1/5)
The optimal bandwidth shrinks as n grows (at rate n^(-1/5)), meaning more data lets you use a finer resolution.
Silverman's normal reference rule (the most widely used starting point):
h_SROT = 0.9 × A × n^(-1/5), where A = min{S, (sample IQR)/1.34}
S is the sample standard deviation
In R: bw = "nrd" invokes Silverman's rule
Other normal reference rules:
h = 1.06 × S × n^(-1/5) (simple normal reference, tends to oversmooth)
h = 0.79 × IQR × n^(-1/5)
More sophisticated methods include cross-validation and plug-in selectors (Sheather-Jones method).
You cannot evaluate f-hat at every possible y, so you choose m plotting points x_1, ..., x_m and connect the values f-hat(x_1), ..., f-hat(x_m) with a smooth curve.
In R, the density() function handles this:
y: data vector
n (or m): number of plotting points (default 512; 100+ is generally enough)
kernel: kernel function ("gaussian" is default)
bw: bandwidth method or numeric value ("nrd" is Silverman's default)
from, to: range of plotting points
Sample size n is very important for density estimation quality
Bandwidth h is the most influential parameter choice
The Gaussian kernel generally produces acceptable results; kernel choice matters far less than bandwidth
The number of plotting points m matters little, provided it is large enough (≥ 100)
When the random variable is strictly positive, the standard KDE can place probability mass below zero. Methods to handle this:
Use from = 0 in R's density() function (simple truncation)
Reflection method: create a mirrored copy of the data at zero, estimate the density on the doubled data, then multiply by 2 for the positive half
Log-KDE: apply the kernel density estimator on the log scale, then transform back (available via the logKDE package in R)
Discrete PMF estimator
f-hat(y_j) = (number of Yi equal to y_j) / n
Histogram density estimator (equal bins)
f-hat(y) = n_j / (n × h) for y in [t_(j-1), t_j)
Histogram density estimator (unequal bins)
f-hat(y) = n_j / (n × h_j) for y in [t_(j-1), t_j)
Kernel density estimator
f-hat(y) = (1/(nh)) × Σ K((y - Yi)/h)
AMISE
AMISE = (1/(nh)) × R(K) + (h⁴/4) × μ_2(K)² × R(f'')
Optimal AMISE bandwidth
h_AMISE = [R(K) / (μ_2(K)² × R(f''))]^(1/5) × n^(-1/5)
Silverman's rule of thumb
h_SROT = 0.9 × min{S, IQR/1.34} × n^(-1/5)
⚠️ Know the difference between a frequency histogram, a relative frequency histogram, and a density histogram. Only the density histogram (relative frequency / class width) is a proper density estimator.
⚠️ With unequal bin widths, plotting raw counts or relative frequencies is misleading. You must divide by the bin width.
⚠️ The KDE formula f-hat(y) = (1/(nh)) × Σ K((y - Yi)/h) is central. Be able to compute individual contributions from each data point (as in the worked example with n = 10, h = 0.3).
⚠️ Bandwidth is the most critical choice in kernel density estimation. Kernel choice is secondary. Know the tradeoff: small h = wiggly (high variance), large h = flat (high bias).
⚠️ Silverman's rule h = 0.9 × A × n^(-1/5) is a common exam formula. Know what A represents.
⚠️ The AMISE has two terms: a variance term proportional to 1/(nh) and a bias-squared term proportional to h⁴. The optimal h balances these.
Q: In a histogram with equal bin widths, what is the difference between plotting frequencies, relative frequencies, and the density estimator?
A: Frequencies are raw counts (n_j). Relative frequencies are proportions (n_j/n). The density estimator divides the relative frequency by the bin width (n_j/(n×h)). All three produce the same shape with equal bins, but only the density estimator has total area = 1.
Q: Why can plotting relative frequencies be misleading with unequal bin widths?
A: Wider bins accumulate more observations simply because they cover a larger range, making those bins appear disproportionately important. Dividing by the bin width corrects for this.
Q: Write the kernel density estimator formula and identify each component.
A: f-hat(y) = (1/(nh)) × Σ K((y - Yi)/h). Here n is the sample size, h is the bandwidth, K is the kernel function, and (y - Yi)/h standardises the distance from y to each data point.
Q: If you use a Gaussian kernel with n = 10 and h = 0.3, what is the contribution of a data point Yi = -0.71 to the estimate at y = -0.8?
A: Compute (y - Yi)/h = (-0.8 - (-0.71))/0.3 = -0.09/0.3 = -0.3. Then (1/(nh)) × K(-0.3) = (1/3) × (1/√(2π)) × exp(-0.045) ≈ 0.127.
Q: What happens to the KDE as bandwidth h increases? As it decreases?
A: As h increases, the estimate becomes smoother (flatter, higher bias, lower variance). As h decreases, the estimate becomes more jagged (lower bias, higher variance, spurious peaks).
Q: State Silverman's rule of thumb for bandwidth selection.
A: h = 0.9 × A × n^(-1/5), where A = min{S, IQR/1.34}, S is the sample standard deviation, and IQR is the interquartile range.
probability density function estimator, PMF estimator, histogram, relative frequency histogram, density histogram, frequency histogram, bin width, class width, number of bins, Scott's rule, Freedman-Diaconis rule, Sturges' rule, kernel density estimator, KDE, bandwidth, smoothing parameter, Gaussian kernel, Epanechnikov kernel, MISE, AMISE, bias-variance tradeoff, Silverman rule, normal reference rule, Sheather-Jones bandwidth, boundary correction, reflection method, logKDE, R density function, ozone data, Principles of Statistics I, Texas A&M, STAT, Tamhane Dunlop