Theoretical Background
Unlike moment-based inequality indicators (e.g., the Gini coefficient), which are highly sensitive to large values in the long tails, indicators which are based solely on quantiles are remarkably resistant to anomalous observations and high distribution skewness.
The core of the inequantiles package is the quantile ratio index (QRI), an indicator that provides a robust measure of inequality, even in small samples, as it considers the entire distribution and is solely based on quantiles.
The Quantile Ratio Index (QRI)
The QRI was introduced by Prendergast and Staudte (2018) as a simple, effective inequality measure of economic inequality. The QRI
- Considers the entire distribution
- Depends solely on quantiles
- Is robust to extreme values
- Does not require a priori choice of specific percentiles
- Is nonparametric
For a continuous non-negative random variable \(Y\) with cumulative distribution function \(F(y)\) and quantile function \(Q(p) = F^{-1}(p)\), \(0 \leq p \leq 1\), define the ratio between symmetric quantiles as:
\[R(p) = \frac{Q(p/2)}{Q(1-p/2)}, \quad 0 \leq p \leq 1.\] \(R(p)\) is a non-decreasing function between 0 and 1, which takes value 1 for any \(0 \leq p \leq 1\) when income (or another economic variable with positive support) is equally distributed among the individuals.
The QRI is then defined as:
\[\text{QRI} = 1 - \int_0^1 R(p) \, dp = 1 - \int_0^1 \frac{Q(p/2)}{Q(1-p/2)} \, dp.\]
The QRI ranges from 0 (perfect equality) to 1 (maximum inequality). It measures the area between the equi-distribution line (\(R(p) = 1\) for all \(p\)) and the actual inequality curve. This can be easily visualised by the following plots:
### Log-Normal distribution
plot_inequality_curve(qfunction = qlnorm, qfun_args = list(meanlog = 9, sdlog = 0.3),
col = "tomato", lty = 2, label = "LogN(sigma=0.3)")
plot_inequality_curve(qfunction = qlnorm, qfun_args = list(meanlog = 9, sdlog = 1.9),
col = "blue", lty = 2, add = TRUE, label = "LogN(sigma=1.9)")
plot_inequality_curve(qfunction = qlnorm, qfun_args = list(meanlog = 9, sdlog = 3.9),
col = "green", lty = 2, add = TRUE, label = "LogN(sigma=3.9)")
For superpopulation models defined by theoretical parametric distributions with known quantile functions, you can compute the exact QRI as:
# Log-normal distribution
superpop_qri(qfunction = qlnorm, meanlog = 9, sdlog = 0.2)
#> [1] 0.2534457
superpop_qri(qfunction = qlnorm, meanlog = 9, sdlog = 1.7)
#> [1] 0.7818321
superpop_qri(qfunction = qlnorm, meanlog = 9, sdlog = 3.2)
#> [1] 0.8781747
# Weibull distribution
superpop_qri(qfunction = qweibull, shape = 2, scale = 30000)
#> [1] 0.522862
superpop_qri(qfunction = qweibull, shape = 1.5, scale = 30000)
#> [1] 0.6003122Consider now a finite population \(U = \{1, \ldots, N\}\), from which a random sample \(s\) of size \(n\) is selected, typically collected with a complex sampling design \(p(s) = Pr(S = s)\), \(\forall s \subseteq U\). Let \(y_j\), \(j \in s\), be the observed values of the variable of interest, with \(y_{(1)}, \ldots, y_{(n)}\) denoting its order statistics. Assume that the sample is drawn according to a certain sampling scheme, with inclusion probability \(\pi_j = Pr(j \in s)\). The corresponding sampling weight \(w_j\) is obtained by the inversion of the inclusion probability, plus, when required, some adjustments for non-response and calibration. Let \(W_j = \sum_{i \in s} w_i \mathbf{1}(i \leq j)\) denote the cumulative sum of weights up to ordered observation \(j\). Let \(\widehat{Q}(p)\) be the \(p\) quantile estimator; how to estimate the quantiles on complex sampling data will be treated in the following sections. For survey data from a finite population, Scarpa, Ferrante, and Sperlich (2025) estimate the QRI using a grid of \(M\) points on \((0, 1)\) as
\[ \widehat{\text{QRI}} = \frac{1}{M} \sum_{m=1}^M \left(1 - \frac{\widehat{Q}(p_m/2)}{\widehat{Q}(1-p_m/2)}\right) \]
where \(p_m = (m-0.5)/M\), for $ m = 1, , M$. By default, \(M = 100\).
\(\widehat{\text{QRI}}\) is strictly sensitive to the choice of the quantile estimator, especially in small samples.
Quantile estimators
The \(p\) quantile estimator can be expressed as a weighted average of order statistics,
\[ \widehat{Q}(p)=y_{(k-1)}+ \left(y_{(k )} - y_{(k- 1)}\right) \left(\frac{p - \widehat{r}_{k - 1 }}{\widehat{r}_{k} - \widehat{r}_{k - 1}} \right), \]
where \(\widehat{r}_{k}\) indicates the estimator of the cdf, namely the plotting position, and the selected order \(k\) is such that \(W_{k-1} - m_{k-1} < W_n p < W_{k} - m_k\), where \(m_k\) is determined by the interpolation method between adjacent data points. Linear interpolation between the points \((\widehat r, y_{(k)})\) gives a quantile estimator for complex sampling data. For \(p=0\) and \(p=1\), define \(\widehat{Q} (0)=y_{(1)}\) and \(\widehat{Q}(1)=y_{(n)}\).
Thecsquantile() function extends standard quantile
estimation through the R function quantile() to survey
data. It adapts the methods described in Hyndman
and Fan (1996) to weighted
data. The possible interpolation rules are summarised in the table below
(see Scarpa, Ferrante, and Sperlich (2025) for further details):
| Estimator | \(\widehat{r}_k\) | \(\widehat{m}_k\) | \(k\) |
|---|---|---|---|
| \(\widehat{Q}_4(p)\) | \(\frac{W_k}{W_n}\) | 0 | \(W_{k-1} \le W_n p \lt W_k\) |
| \(\widehat{Q}_5(p)\) | \(\frac{W_k-\frac{1}{2}w_k}{W_n}\) | \(\frac{w_k}{2}\) | \(W_{k-1} - \frac{w_{k-1}}{2} \le W_n p \lt W_{k} - \frac{w_{k}}{2}\) |
| \(\widehat{Q}_6(p)\) | \(\frac{W_k}{W_n+w_n}\) | \(w_np\) | \(W_{k-1} \le (W_n + w_n)p \lt W_{k}\) |
| \(\widehat{Q}_7(p)\) | \(\frac{W_{k-1}}{W_{n-1}}\) | \(w_k - w_np\) | \(W_{k-2} \le W_{n-1}p \lt W_{k-1}\) |
| \(\widehat{Q}_8(p)\) | \(\frac{W_k-\frac{1}{3}w_k}{W_n+rac{w_n}{3}}\) | \(\frac{w_k}{3} + \frac{w_n}{3}p\) | \(W_{k-1} - \frac{w_{k-1}}{3} \le (W_{n} - \frac{w_n}{3})p \lt W_{k} - \frac{w_k}{3}\) |
| \(\widehat{Q}_9(p)\) | \(\frac{W_k-\frac{3}{8}w_k}{W_n+rac{1}{4}w_n}\) | \(\frac{3}{8}w_k + \frac{w_n}{4}p\) | \(W_{k-1} - \frac{3w_{k-1}}{8} \le (W_{n} + \frac{w_{n}}{4})p \lt W_{k} - \frac{3w_{k}}{8}\) |
They can be easily computed on survey data, as:
# Compute weighted quartiles
csquantile(y = synthouse$eq_income,
weights = synthouse$weight,
probs = c(0.25, 0.5, 0.75),
type = 4)
#> 25% 50% 75%
#> 12353.29 20014.13 32222.45
csquantile(y = synthouse$eq_income,
weights = synthouse$weight,
probs = c(0.25, 0.5, 0.75),
type = 7)
#> 25% 50% 75%
#> 12353.29 20019.77 32224.95
# Compare without considering the sampling weights
csquantile(synthouse$eq_income, probs = c(0.25, 0.5, 0.75), type = 4)
#> 25% 50% 75%
#> 12910.30 20428.45 32528.88An extension of the Harrell-Davis estimator to survey data is also provided, as \(\widehat{Q}_{HD}(p)=\sum_{j \in s} \widehat{\mathcal{W}}_{j}(p) y_{(j)}\), where
\[ \begin{aligned} \widehat{\mathcal{W}}_{j}(p) = & \, b_{(W_{j} / W_n)}\{(W_n+ w_n) p, W_n - (W_n+ w_n)p + w_n\} \\ & - b_{(W_{j - 1}/ W_n)}\{(W_n+ w_n) p, W_n - (W_n+ w_n)p + w_n\} \end{aligned} \]
# Harrell-Davis weighted median
csquantile(y = synthouse$eq_income,
weights = synthouse$weight,
probs = 0.5,
type = "HD")
#> 50%
#> 20017.37Differences among the quantile estimators are particularly evident in small samples and in the distribution tails, as demonstrated in the following example:
# Compare different quantile types by NUTS3
types <- c(4, 5, 6, 7, 8, 9, "HD")
areas <- unique(synthouse$NUTS3)
# Function to compute QRI for all types in one area
compare_quantiles <- function(region_code, data = synthouse) {
idx <- which(data$NUTS3 == region_code)
results <- sapply(types, function(t) {
csquantile(y = data$eq_income[idx],
weights = data$weight[idx],
type = t,
probs = 0.95)
})
return(results)
}
# Compute for all areas
results_quantiles <- sapply(areas, compare_quantiles)
rownames(results_quantiles) <- types
colnames(results_quantiles) <- areas
## === Quantile estimators for Each NUTS3 Region ==="
print(head(t(results_quantiles), n = 10))
#> 4 5 6 7 8 9 HD
#> N01005 70988.83 75538.32 87604.07 69761.89 80343.39 79093.40 76388.02
#> NE06004 50988.47 50988.47 50988.47 50992.05 50988.47 50988.47 50988.47
#> N05003 81116.08 82948.72 81211.01 83560.96 82496.53 82633.53 81370.54
#> NE06002 47742.34 47749.81 47749.81 47749.81 47749.81 47749.81 47749.81
#> N02004 52599.79 53137.63 52644.75 53559.00 52883.63 52956.29 52818.39
#> NO02004 66674.78 66674.78 67773.72 61032.04 66904.21 66771.79 66674.78
#> NO03002 67724.83 71393.13 71459.35 71351.81 71411.51 71406.64 72491.41
#> N06004 45313.69 45317.87 45313.69 45706.20 45313.69 45313.69 45313.69
#> NE02004 57268.49 57299.09 57303.00 57297.80 57299.86 57299.64 57279.50
#> S04005 58411.94 62602.84 60288.65 62670.74 61561.74 61786.33 62670.69The package supports multiple quantile estimation methods (types 4-9
and Harrell-Davis) into quantile-based inequality indicators estimators.
Rule type = 6is recommended for the QRI, as it is been
shown by Scarpa, Ferrante, and Sperlich (2025) to lead to the least
biased estimates.
QRI estimator sampling variance
For complex surveys, the rescaled bootstrap method (Rao and Wu 1988; Rao, Wu, and Yue 1992) is recommended for the estimation of the QRI estimator sampling variance, as demonstrated by Scarpa, Ferrante, and Sperlich (2025). This can be easily estimated on data collected through two-stage stratified sampling design as
Other Inequality Indicators
Most common quantile-based inequality indicators in the literature
are the quintile share ratio (QSR), the Palma ratio and the interdecile
ratios. All estimators implemented in this package allow the user to
choose the quantile estimation method via the type
argument.
Palma Ratio
The Palma ratio compares the top 10% to the bottom 40% aggregated income: \[ \widehat{{Palma}} = \frac{\sum_{j \in s}w_j y_j \mathbf{1}\left\{ y_j \geq \widehat{Q}(0.9)\right\} }{\sum_{j \in s} w_j y_j\mathbf{1}\left\{ y_j \leq \widehat{Q}(0.4)\right\} } \ . \]
Percentile Ratios
Very often National Statistical Offices measure inequality using ratios between percentiles, for example: \[ \widehat{{P}90/{P}10} = \frac{\widehat{Q}(p=0.9)}{\widehat{Q}(p=0.1)} \]
Default the P90/P10 is estimated.
Comparing Indicators Across Groups
Multiple inequality indicators exist in the literature, each capturing different aspects of the economic distribution. Different indicators may provide complementary perspectives on inequality, so researchers should choose the most appropriate measure based on their specific research questions and the distributional features of interest.
The inequantiles() function computes all quantile-based
indicators at once, or any subset of them, among QRI, QSR, Palma ratio
and interquantiles ratio.
# Compare all indicators
inequantiles(y = synthouse$eq_income,
weights = synthouse$weight,
indicators = "all")
#> Quantile-based inequality indicators
#> -------------------------------------
#> Estimate
#> qri 0.5691
#> qsr 7.0161
#> palma 1.5787
#> P90/P10 5.9206
#> gini 0.3715Standard errors can be obtained via the rescaled bootstrap by setting
se = TRUE. Note that this requires specifying the survey
design variables (strata, and optionally psu
and N_h).
# QRI and QSR with their standard errors for one region via rescaled bootstrap
nord <- synthouse[synthouse$NUTS1 == "N", ]
qri_qsr <- inequantiles(
y = nord$eq_income,
weights = nord$weight,
indicators = c("qri", "qsr"),
se = TRUE,
data = nord,
strata = "NUTS2",
psu = "municipality",
B = 200,
seed = 42
)