Sobol' indices

For random or stochastic models, one may wonder which input provides more information to the model. We will assume that our model is of the form $\color[rgb]{1,1,.99} Y=f(X_1,\dots,X_d)$, where $\color[rgb]{1,1,.99} Y$ is our random variable, a function of $\color[rgb]{1,1,.99} d$ independent random variables. From the integration point of view, we can answer these questions to measure the effective dimensionality of our function and understand why quasi-Monte Carlo cubatures might perform really well/badly for some high dimensional functions.

Sobol' indices are one measure that can answer this question. To define these indices, we need to first know the definition of the ANOVA decomposotion. For $\color[rgb]{1,1,.99} f\in L^2\left([0,1]^d\right)$, and $\color[rgb]{1,1,.99}\mathcal{D}=\{1,\dots,d\}$,

{\color[rgb]{1,1,.99} \begin{align*} f(\vec{x})=\sum_{u \subseteq \mathcal{D}} f_u(\vec{x})\, , \qquad f_{\varnothing} = \mu\, , \end{align*} }
where
{\color[rgb]{1,1,.99} \begin{align*} f_u(\vec{x})= \int_{[0,1]^{d-|u|}} f(\vec{x}) d{\vec{x}}_{-u} - \sum_{v \subset u} f_v(\vec{x})\, . \end{align*} }

In the notation above, $\color[rgb]{1,1,.99} |u|$ is the cardinality of $\color[rgb]{1,1,.99} u$, and $\color[rgb]{1,1,.99} -u:=u^c=\mathcal{D}\setminus u$.

If we define the variances to be:

{\color[rgb]{1,1,.99} \begin{align*} \sigma^2_{\varnothing} = 0\, ,\qquad \sigma^2_{u} =\int_{[0,1]^{d}} f_u(\vec{x})^2 d{\vec{x}}\, ,\qquad \sigma^2 =\int_{[0,1]^{d}} \left(f(\vec{x})-\mu\right)^2 d{\vec{x}}\, . \end{align*} }

then, one can easily obtain the ANOVA identity,

{\color[rgb]{1,1,.99} \begin{align*} \sigma^2 = \sum \limits_{u \subseteq\mathcal{D}} \sigma_u^2 \, . \end{align*} }

In 91, Sobol' introduced the global sensitivity indices,

{\color[rgb]{1,1,.99} \begin{align*} \underline{\tau}_u^2 = \sum_{\substack{v \subseteq u \\ v\,\in\mathcal{D}}} \sigma_v^2\, , \quad \text{ and } \quad \overline{\tau}_u^2 = \sum_{\substack{v \cap u\neq\varnothing \\ v\,\in\mathcal{D}}} \sigma_v^2\, . \end{align*} }

Using plain words, $\color[rgb]{1,1,.99} \underline{\tau}_u^2$ measures the variance explained only by input dimensions $\color[rgb]{1,1,.99} u$, while $\color[rgb]{1,1,.99} \overline{\tau}_u^2$ measures the variance explained by $\color[rgb]{1,1,.99} u$ and its interactions with any other set of dimensions. It may be easier to understand this definition considering the properties $\color[rgb]{1,1,.99} \underline{\tau}_u^2\leq \overline{\tau}_u^2$ and $\color[rgb]{1,1,.99} \underline{\tau}_u^2 + \overline{\tau}_{-u}^2 =\sigma^2$.

Some people also work with the normalized version of these indices which is $\color[rgb]{1,1,.99} S_u:=\underline{\tau}_u^2/\sigma^2$ and $\color[rgb]{1,1,.99} \overline{S}_u:=\overline{\tau}_u^2/\sigma^2$. From the probabilistic point of view, considering each input variable to be independent and uniformly distributed, Sobol' indices are also written under this equivalent form (with $\color[rgb]{1,1,.99} u$ indicating the dimensions that are not frozen):

{\color[rgb]{1,1,.99} \begin{align*} S_{u}=\frac{{\rm Var} \left[ \mathbb{E} \left(f({\vec{X}})|{\vec{X}}_{u}\right)\right] }{{\rm Var}\left(f({\vec{X}})\right)}= 1-\frac{\mathbb{E}\left[ {\rm Var} \left(f({\vec{X}})|\vec{X}_u\right)\right]}{{\rm Var}\left(f({\vec{X}})\right)}, \end{align*} }
{\color[rgb]{1,1,.99} \begin{align*} \overline{S}_{u}=1-\frac{{\rm Var} \left[ \mathbb{E} \left(f({\vec{X}})|\vec{X}_{-u}\right)\right]}{{\rm Var}\left(f({\vec{X}})\right)}=\frac{\mathbb{E}\left[ {\rm Var} \left(f({\vec{X}})|{\vec{X}}_{-u}\right)\right] }{{\rm Var}\left(f({\vec{X}})\right)}. \end{align*} }

Because these indices are functions over integrals, we can carefully define new estimators to adapt our automatic cubatures to estimate them (ArXiv:1411.1966, ArXiv:1410.8615).

Simple example

The Brownian motion is widely used in many applications. Thus, the following example might be interesting to understand what Sobol' indices can explain in some of these applications. For instance, if we want to estimate the expected value of the maximum of a discretized Brownian motion, (according to our notation in BM example, and $\color[rgb]{1,1,.99} \Sigma$ the covariance matrix of the BM), our goal is to estimate

{\color[rgb]{1,1,.99} \begin{align*} \mathbb{E} \left(\max(\vec{B}_t)\right) &= \int_{\mathbb{R}^d} \max(\vec{B}_t) \frac{{\rm e}^{\vec{B}_t^T\Sigma^{-1}\vec{B}_t}}{(2\pi)^{d/2}|\Sigma|^{1/2}}\,d\vec{B}_t \\ &= \int_{[0,1]^d} \max(f_{\rm Chol}\vec{x}) \,d\vec{x} \\ &= \int_{[0,1]^d} \max(f_{\rm PCA}\vec{x}) \,d\vec{x}. \end{align*} }

The last two equalities are obtained through two different subsitutions corresponding to the Cholesky construction of the Brownian motion, and the PCA construction. For instance, if we take 10 dimensions and $\color[rgb]{1,1,.99} t_i=i/10$, using the Cholesky construction $\color[rgb]{1,1,.99} S_{1} = 24\%$, $\color[rgb]{1,1,.99} S_{2} = 16\%$, $\color[rgb]{1,1,.99} S_{3} = 13\%$, $\color[rgb]{1,1,.99} S_{4} = 10\%$, and $\color[rgb]{1,1,.99} \overline{S}_{1}=24\%$, $\color[rgb]{1,1,.99} \overline{S}_{2}=18\%$, $\color[rgb]{1,1,.99} \overline{S}_{3}=15\%$, $\color[rgb]{1,1,.99} \overline{S}_{4} = 13\%$. However, using the PCA construction, $\color[rgb]{1,1,.99} S_{1} = 89\%$, $\color[rgb]{1,1,.99} S_{2} = 2\%$, $\color[rgb]{1,1,.99} S_{3} = 1\%$, $\color[rgb]{1,1,.99} S_{4} = 0\%$, and $\color[rgb]{1,1,.99} \overline{S}_{1}=94\%$, $\color[rgb]{1,1,.99} \overline{S}_{2}=7\%$, $\color[rgb]{1,1,.99} \overline{S}_{3}=3\%$, $\color[rgb]{1,1,.99} \overline{S}_{4} = 2\%$. This shows that the PCA construction will be better to estimate the above expectation using quasi-Monte Carlo methods since it has a lower effective dimensionality.

For those who were curious about the actual value of $\color[rgb]{1,1,.99}\mathbb{E} \left(\max(\vec{B}_t)\right)$, it is approximately 0.5935. Nevetheless, the Brownian Motion is a continuous time process. Thus, in order to estimate the expected value on a continuous time Browninan motion, one could either use a multilevel method, or a multivariate decomposition method to estimate this infinite dimensional integral.