A sensitivity analysis of the PAWN sensitivity index

TL;DR

This paper investigates how the PAWN sensitivity index's robustness is affected by its design parameters, revealing potential biases and challenges in reliably ranking model inputs compared to Sobol' indices.

Contribution

It provides a comprehensive sensitivity analysis of PAWN's design parameters and compares its robustness with Sobol' indices, highlighting limitations in input influence detection.

Findings

PAWN's results are sensitive to its design parameters.

Design uncertainties can bias PAWN's input ranking.

PAWN may struggle to distinguish influential from non-influential inputs.

Abstract

The PAWN index is gaining traction among the modelling community as a sensitivity measure. However, the robustness to its design parameters has not yet been scrutinized: the size ($N$) and sampling ($\varepsilon$) of the model output, the number of conditioning intervals ($n$) or the summary statistic ($\theta$). Here we fill this gap by running a sensitivity analysis of a PAWN-based sensitivity analysis. We compare the results with the design uncertainties of the Sobol' total-order index ($S_{Ti}^*$). Unlike in $S_{Ti}^*$, the design uncertainties in PAWN create non-negligible chances of producing biased results when ranking or screening inputs. The dependence of PAWN upon ($N,n,\varepsilon, \theta$) is difficult to tame, as these parameters interact with one another. Even in an ideal setting in which the optimum choice for ($N,n,\varepsilon, \theta$) is known in advance, PAWN might not allow to distinguish an influential, non-additive model input from a truly non-influential model input.

Tables

Table 1. Table 1: Summary of the parameters and their distribution for both PAWN ( N , n , ε , θ 𝑁 𝑛 𝜀 𝜃 N,n,\varepsilon,\theta ) and S T ​ i ∗ superscript subscript 𝑆 𝑇 𝑖 S_{Ti}^{*} ( N , θ 𝑁 𝜃 N,\theta ). D ​ U 𝐷 𝑈 DU stands for discrete uniform.

Parameter	Description	Distribution
$N$	Total number of runs	$𝒰(200,2000)$
$n$	Number of conditioning intervals	$𝒰(5,20)$
$\epsilon$	Randomness in the sampling of $Y_U$	$𝒰(1,{10}^3)$
$\theta$	Summary statistic/estimator	$𝒟𝒰(1,3)$

Equations

$$Y=X_1/X_2$$

$$Y=\sin (X_1)+a\sin (X_2{)}^2+b{X}_3^4\sin (X_1)$$

$$Y=\prod _{i=1}^k\frac{|4X_i-2|+a_i}{1+a_i}$$

$$Y=$$

$$+\sum _{i<j<l}^{20}{\beta }_{i,j,l}{w}_i{w}_j{w}_l$$

$$+\sum _{i<j<l<s}^{20}{\beta }_{i,j,l,s}{w}_i{w}_j{w}_l{w}_s$$

$$V_{i}=V_{X_{i}}\big{[}E_{\textbf{X}_{\sim i}}(Y|X_{i})\big{]}$$

$$S_i=\frac{V_i}{V_Y}$$

$$V_Y=\sum _{i=1}^k{V}_i+\sum _i\sum _{i<j}{V}_{ij}+...+V_{1,2,...,k}$$

$$V_{ij}=$$

$$\displaystyle-V_{X_{i}}\big{[}E_{\textbf{X}_{\sim i}}(Y|X_{i})\big{]}$$

$$\displaystyle-V_{X_{j}}\big{[}E_{\textbf{X}_{\sim j}}(Y|X_{j})\big{]}$$

$$S_{ij}=\frac{V_{ij}}{V_Y}$$

$$S_{Ti}=\frac{E_{X_{\sim i}}\big{[}V_{X_{i}}(Y|X_{\sim i})\big{]}}{V_{Y}}$$

Full text

Keywords: Uncertainty, Environmental Modelling, Statistics, Risk

\saythanks

1 Introduction

[1, 2] have recently published in Environmental Modelling & Software a new measure for sensitivity analysis, the PAWN index. Like other moment-independent approaches (i.e. entropy-based [3], density-based [4, 5]), PAWN does not resort to statistical second-order moments such as variance to apportion output uncertainty to the model parameters. Instead, it relies on Cumulative Distribution Functions (CDFs) to characterize the maximum distance between the unconditional output distribution $Y_{U}$, i.e. obtained by moving all parameters simultaneously, and the conditional output distribution $Y_{C_{ij}}$, i.e. obtained by fixing the $i$-th parameter to $j=1,2,...,n$ values or intervals within its uncertainty range. The difference between $Y_{U}$ and $Y_{C_{ij}}$ is assessed via the Kolmogorov-Smirnov test, although other distance-based tests, such as the Anderson-Darling’s, may also be used [6]. The final PAWN index for a given parameter is obtained by calculating the mean, the median, the maximum or any other summary statistic over all the KS values computed between $Y_{U}$ and $Y_{C_{ij}}$.

The most up-to-date approximation to the PAWN index, named \saythe generic approach [2], is as follows: let there be an A matrix with $v=1,2,...,N$ rows and $i=1,2,...,k$ parameters. After computing the model output $Y$, the range of variation of the $i$-th parameter is split into $j=1,2,...,n$ intervals of size $N_{c}$ (where $N_{c}\approx N/n$). The model output linked to the $j$-th interval is used as the conditional model output $Y_{C_{ij}}$. The unconditional model output $Y_{U}$ can concur with the whole model output or can be a random sub-sample of the same size as $N_{c}$. With this approach, the total number of model runs to compute PAWN is fully determined by $N$ [2].

Based on trials with the [3] function, the [7] function, the SWAT model [6] or a wind-energy converter model [8], it has been observed that PAWN might reach convergence much faster than Sobol’ indices. A key question, however, is to know how the selection of $N,n$, the sampling of $Y_{U}$ or the summary statistic affects the accuracy of PAWN. This is timely given the widespread adoption of the index: since its inception in 2015, PAWN has been cited 94 times, with the number of citations stably increasing from 5 in 2015 to 31 in 2019. Most of the works quoting PAWN are from the environmental sciences (53), followed by engineering (24) and computer science (23) (Scopus search on September 25 2019). Gaining a systematic insight into the internal functioning of PAWN shall thus allow the modelling community to better appraise its robustness, thus increasing its transparency as well as our awareness of its advantages and limitations.

Here we assess the sensitivity of PAWN to the main structural uncertainties involving its calculation, an exercise that might be termed \saya sensitivity analysis of a sensitivity analysis (SA of SA). This expression was used by [9] to study how sensitivity indices are affected by uncertainties in the probability distributions used to describe the model inputs, a work that actually falls into the tradition of \sayprobability of probability or \sayprobability of frequency described in [10] (see also [11]). A similar analysis of the sensitivity of results to changes in the range of input factors can be found in [12]. Other approaches include exploring several sensitivity measures [13, 14, 15]. More recently, [16] included considerations of form of the output, sample size, choice of method, measuring interactions, range of distributions, inclusion or exclusion of non-behavioural runs (simulations), and use of dummy variables.

Our experiment departs from previous studies in that it assesses the uncertainties embedded in the structural design of sensitivity indices, PAWN in that case. We thus explore the implications of the last ring in the chain of uncertainties characterizing any sensitivity analysis, the last stretches of the garden of forking paths that need to be crossed to arrive at a result in any modelling exercise [17, 18]. Exploring all sources of uncertainty is indeed an important pursuit, including the notions of context, purpose and motivations as suggested in sensitivity auditing [19]. However, we believe that the papers just reviewed [9, 15, 12, 16] are not a SA of a SA but more general instances of uncertainty exploration. We take the meaning of SA of SA literally (i.e. the analysis of the sensitivity of a sensitivity analysis to its own design parameters), and thus exclude in this paper the exploration of sensitivities to changing methods, assumptions or model designs, however worthy these analyses might be.

2 Materials and methods

In order to cross-check our approach, we match the design uncertainties of PAWN against the design uncertainties of the Sobol’ total-order index ($S_{Ti}^{*}$), a well-established measure to determine how much a given model parameter interacts with the rest111The asterisk in $S_{Ti}^{*}$ is used to differentiate between Sobol’ total-order index ($S_{Ti}^{*}$), i.e. the index against which PAWN is matched, and the index used to assess the extent to which $S_{Ti}^{*}$ interacts with its design parameters ($N,\theta$), $S_{Ti}$.. Although different in scope and nature, the fact that unconditional CDFs are also affected by interactions paves the way for the uncertainties in PAWN and $S_{Ti}^{*}$ to be explored in parallel.

For PAWN, we focused on four uncertain parameters: the total number of model runs ($N$), the number of conditioning intervals ($n$), the randomness derived from the sampling of the unconditional model output $Y_{U}$ ($\varepsilon$), and the summary statistic ($\theta$). For $S_{Ti}^{*}$, we focused on two: $N$ and $\theta$, the latter reflecting the different estimators existing to compute the total-order effect [20, 21, 22] (Table 1). We assessed how different combinations of values for these uncertain parameters condition the PAWN/$S_{Ti}^{*}$ index using four different test functions that yield a skewed model output (Fig. 1): the [3]’s, which reads as

[TABLE]

where $X_{1}\sim\chi^{2}(10)$ and $X_{2}\sim\chi^{2}(13.978)$; the [7]’s, which reads as

[TABLE]

where $a=2,b=1$ and $(X_{1},X_{2},X_{3})\sim\mathcal{U}(-\pi,+\pi)$; the Sobol’ G [23]’s, which reads as

[TABLE]

where $k=8$, $X_{i}\sim\mathcal{U}(0,1)$ and $a=(0,1,4.5,9,99,99,99,99)$; and the [24] function, which reads as

[TABLE]

where $w_{i}=2(X_{i}-0.5)$ for all $i$ except for $i=4,5,7$, where $w_{i}=2(1.1X_{i}/(X_{i}+0.1)-0.5)$, $X_{i}\sim\mathcal{U}(0,1)$, and

$\beta_{i}=20,i=1,2,...,10$,

$\beta_{i,j}=-15,i=1,2,...,6$,

$\beta_{i,j,l}=-10,i=1,2,...,5$,

$\beta_{i,j,l,s}=5,i=1,2,...,4$ [25].

We chose these four skewed benchmark functions in order to provide PAWN with a favourable test ground: density-based measures might better characterize the uncertainty in skewed model outputs than variance-based measures as the former do not rely on any statistical moment [26].

We selected the distributions of $(N,n,\varepsilon,\theta)$ based on previous work on PAWN and the $S_{Ti}^{*}$ index [1, 2, 27] and some preliminary tests. Firstly, we observed that $N\approx 2000$ was sufficient for the PAWN index of many model inputs to be very close to convergence. Sobol’ $S_{Ti}^{*}$, on the other hand, required a larger number of runs (Fig. S1). This corroborated previous observations on the faster convergence rate of PAWN compared to Sobol’ indices [1, 6]. By defining $N\sim\mathcal{U}(200,2000)$ we set our study in a scenario where the uncertainty with regards to the required sample size needed to compute robust PAWN indices is moderate. Indeed, this is usually the case: the analyst might have a cap on the total number of model runs available with the computing resources at hand, but no prior information on the minimum sample size required to ensure stable sensitivity indices for all model inputs.

We defined the distribution of $n$ based on [2], who suggested to start with $n=10$ and vary $n$ some units up and down to check its effect (provided that $n>5$).

For $\varepsilon$, we set $10^{3}$ different starting points (seeds) for the pseudo-random number sequence used to generate the indices (from 1 to $N$) to sample $Y_{U}$. This ensured 1) negligible chances of the same seed overlapping with the same value for $N$, thus introducing determinism into a process that should be mainly stochastic, and 2) that the randomness in the sampling of $Y_{U}$ is assessed in terms of its relative influence in the computation of PAWN.

With regards to $\theta$, for PAWN we used the mean, the median and the maximum as a summary statistic for the KS values when $\theta=1,\theta=2$ and $\theta=3$. For $S_{Ti}^{*}$, we used the estimators by [20], [21] and [22] when $\theta=1,\theta=2$ and $\theta=3$ respectively.

In order to estimate the uncertainty propagated by ($N,n,\varepsilon,\theta$) (resp. $N,\theta$) to PAWN (resp. $S_{Ti}^{*}$), we created a $(2^{13},2k)$ sample matrix for each function using Sobol’ quasi-random number sequences, where $k=4$ (resp. $k=2$), and transformed the columns into their appropriate distributions (Table 1). The first $k$ matrix was labelled A and the second $k$ matrix, B. Our model ran row-wise in both the A and B matrices, as follows: based on the information contained in the $v$-th row, it created a Sobol’ matrix of size $N^{(v)}$ for PAWN [of size $int(N^{(v)}/(k+1))$ for $S_{Ti}^{*}$], and computed either the [3], the [7], the Sobol’ G [23] or the [24] function. Then, for each model input $i$, it either calculated the PAWN*(v)* index following the conditions set by $n^{(v)}$, $\varepsilon^{(v)}$ and $\theta^{(v)}$, or the $S_{Ti}^{*(v)}$ index according to $\theta^{(v)}$.

We estimated how sensitive PAWN ($S_{Ti}^{*}$) indices are to uncertainty in ($N,n,\varepsilon,\theta$) (resp. $N,\theta$) by means of Sobol’ indices [28]. For a model of the form $Y=f(X_{1},X_{2},...,X_{k})$, where $Y$ is a scalar and $X_{1},X_{2},...,X_{k}$ are independent parameters described by known probability distributions, we can measure how sensible $Y$ is to a given parameter $X_{i}$ with

[TABLE]

where $E_{\textbf{X}_{\sim i}}(Y|X_{i})$ is the expected value of $Y$ calculated over all possible values of all parameters except the $i$-th, which is kept fixed. By dividing Equation 5 by the unconditional model output variance, we obtain the first order sensitivity index for $X_{i}$, which describes the proportion of variance in the model output caused by $X_{i}$:

[TABLE]

We can then decompose the unconditional model output variance $Y$ as the sum of conditional variances up to the $k$-th order:

[TABLE]

where

[TABLE]

From this, we can derive the second-order index $S_{ij}$, which explains the proportion of variance due to the interaction between $X_{i}$ and $X_{j}$:

[TABLE]

and so on until order $k$. However, estimating all terms in Equation 7 is unattainable when $k$ is large, as they result in $2^{k}-1$. In this case, we can compute the total order index or $S_{Ti}$, which measures the proportion of variance due to the first-order effect of $X_{i}$ jointly with its interactions with the other parameters [21]:

[TABLE]

For PAWN, since $k=4$, we computed first ($S_{i}$), second ($S_{ij}$), third ($S_{ijk}$) and total-order ($S_{Ti}$) Sobol’ indices of ($N,n,\varepsilon,\theta$). For $S_{Ti}^{*}$, since $k=2$, we just computed $S_{i}$ and $S_{Ti}$. In both settings we used the [27] and the [20] estimators to compute $S_{i}$ and $S_{Ti}$ respectively, as per the established best practices. All the workflow is summarised in Fig. S2 and the $R$ code to replicate our results is available in GitHub.

3 Results

3.1 Uncertainty analysis

Fig. 2 presents the uncertainty distribution of PAWN and $S_{Ti}^{*}$ for each model input and function. The results can be matched against Fig. S3, where we display how PAWN and Sobol’ indices look like once their design parameters are fixed and the total number of model runs is set at $N=4000$.

For PAWN, we ran three simulations (Fig. 2a):

With $\theta$ including the mean, the median and the maximum as possible summary statistics ($max\in\theta$ setting). 2. 2.

With $\theta$ including the mean and the median only ($max\notin\theta$ setting). This aimed at isolating the effect that extreme KS values, which might be obtained for specific conditioning intervals, have in the final PAWN index. 3. 3.

With $N\sim\mathcal{U}(2500,4000)$, $n\sim\mathcal{U}(15,20)$ and ($max\notin\theta$) (\sayOptimum setting). This latter run reflects an ideal scenario, one in which the number of model runs needed to achieve convergence is known in advance and the uncertainty in $n$ and $\theta$ is reduced to the minimum expression.

For $S_{Ti}^{*}$ we run one simulation to assess its sensitivity to the uncertainty in the total number of model runs and the estimator ($N,\theta$) [with $N\sim\mathcal{(}200,2000)$] (Fig. 2b). Sobol’ indices might eventually take on values outside the range $[0,1]$ due to numerical artifacts created during the computation. This is widely known among sensitivity analysts and managed by considering $(S_{Ti}^{*}<0)\approx 0$ and $(S_{Ti}^{*}>1)\approx 1$. Fig. 2b only presents $S_{Ti}^{*}\in[0,1]$ to allow for a better comparison with the values produced by the PAWN index. A plot showing the distribution of values with $S_{Ti}^{*}\notin[0,1]$ is presented in Fig. S4.

3.1.1 Factor prioritization

Fig. 2 shows that, in a factor prioritization context, i.e. when the aim is to sort the parameters according to their contribution to the model output variance [29], the uncertainty in the value of the design parameters might cause model inputs to overlap, thus raising the likelihood of producing a biased ranking. In order to get precise figures for this overlap, we computed the coefficient of overlapping, i.e. the area lying under the density curves of two different model inputs, following [30] (see Fig. S5 for a presentation of Fig 2 with density curves instead).

In the case of the [7] function, the overlap between $X_{2}$ and $X_{3}$ is of 10% if PAWN is used under $max\in\theta$ (4% under $max\notin\theta$), despite $X_{2}$ being non-influential. The overlap with $S_{Ti}^{*}$ is 0. With the [23] G function, PAWN mistakes $X_{1}$ for $X_{2}$ and $X_{2}$ for $X_{3}$ 11% and 20% of the time respectively if $max\in\theta$ (4% and 10% if $max\notin\theta$), and might even bias the ranking of $X_{2}$ and $X_{4}$ (14% and 9% overlap in $max\in\theta$ and $max\notin\theta$ respectively). On the contrary, if $S_{Ti}^{*}$ is used, the overlap between $X_{1}-X_{2}$, $X_{2}-X_{3}$ and $X_{2}-X_{4}$ is $\sim$0%, 2% and 1% respectively.

The case of the [3] function deserves a specific comment. Fig. 2 shows that 29% of $X_{1}$ and $X_{2}$ values overlap if using PAWN under $max\in\theta$ (13% if $max\notin\theta$), whereas for $S_{Ti}^{*}$ the degree of overlap is 55%. [3], however, stated that $X_{1}$ was more influential than $X_{2}$ based on relative-entropy sensitivity methods. [1] used the [3] function to back up their claim of PAWN outperforming Sobol’ indices due to the former being able to discriminate the higher influence of $X_{1}$ much better than $S_{Ti}^{*}$. However, the analytical values of $S_{Ti}^{*}$ in the [3] function are identical at $X_{1}=X_{2}=0.546$ (see Table 2 in [3]). This means that the overlap between $X_{1}$ and $X_{2}$ shown in Fig. 2b is to be expected and does not result from the sensitivity of $S_{Ti}^{*}$ to the uncertainties in its own structural design.

3.1.2 Factor screening

In a factor screening context, i.e. when the aim is to distinguish influential from non-influential parameters [29], the uncertainty in the PAWN design parameters might also lead to erroneous results, regardless of whether $max\in\theta$ or $max\notin\theta$. In the case of the [23] G function, 10-14% (7-8%) of the probability density of $X_{2}$, which has a non-nihil effect, overlaps with $X_{4},...,X_{8}$, which have no effect at all. The overlap between influential and non-influential model inputs if $S_{Ti}^{*}$ is used, on the other side, ranges between $\sim$0-2%.

The poor screening power of PAWN stands out in the case of the [24] function: 18-25% (1-8%) of the probability density of $X_{8},...,X_{10}$, which have a moderate first-order effect, overlap with that of $X_{11},...,X_{20}$, which are non-influential (Fig. S2A). In the case of $S_{Ti}^{*}$, the overlap ranges between 0-5%. The chances of PAWN mistaking relevant for non-relevant parameters is even higher in the case of parameters whose influence in the model output is through interactions only: this is the case of $X_{1},...,X_{6}$ in the [24] function, whose degree of overlap with $X_{11},...,X_{20}$, which are non-influential, range between 75-95% (Figs. 2–S3). In all these cases, the uncertainty in the design parameters of $S_{Ti}^{*}$ leads to a 0-4% overlap between influential and non-influential parameters.

Fig. 2 also shows that in the \sayOptimum setting the overlap between the model inputs is considerably reduced for PAWN. In the case of the [3] and the [7] functions, the percentage of overlap goes down to zero. However, the chances of wrongly ranking/screening the model inputs remain non-negligible for both the [23] G and the [24] functions. In the former, there is 17–28% overlap between the slightly influential model input $X_{3}$ and $X_{4},...,X_{8}$, whose effect can not be differentiated from the approximation error. With regards to the [24] function, the chances of characterizing as non-influential parameters that have a significant non-additive effect in the model output remain very high: the overlap of $X_{1},...,X_{6}$ with $X_{11},...,X_{20}$, for instance, range between 35–90% (Figs. 2). As shown in Fig S3, the volatility in the computation of PAWN does not allow to distinguish $X_{11},...,X_{20}$ from a dummy, non-influential model input.

3.2 Sensitivity analysis

Figs. 3 presents the Sobol’ first ($S_{i}$) and total ($S_{Ti}$) indices for each of the settings of our analysis after pooling the values from all functions and parameters [the Sobol’ indices for each function and design parameter are shown in Figs. S6–S7 (PAWN), and Figs. S8–S9 ($S_{Ti}^{*}$)]. Figs. 3 thus inform on how much each design parameter contribute uncertainty to each sensitivity index. In Fig. S10 of the Supplementary Information file, we prove that the results displayed in Fig. 3 are robust without having to offset the stronger weight that the Morris function might have in defining the trends due to its much larger number of model inputs.

As shown in Fig. 3a, the first-order effect of $\theta$ and $N$ is much more variable in the $max\in\theta$ setting, suggesting that their degree of contribution to the PAWN index uncertainty might considerably be function-dependent. This variability is highly reduced in the $max\notin\theta$ setting, which provides a more robust account of the extent to which each design parameter contributes to define the PAWN index. In this setting, the selection of the initial sample size ($N$) and the number of conditioning intervals ($n$) convey up to c. 60% and c. 20% of the PAWN index uncertainty respectively. The stochasticity in the sampling of $Y_{U}$ ($\varepsilon$) is mostly influential through interactions in both settings, whereas the selection of the summary statistic ($\theta$) has a nearly nihil effect in the $max\notin\theta$ setting. Remarkably, interactions are also significant in the \sayOptimum setting, especially those involving $\varepsilon$, $n$ and $N$.

To gain further insights into the structure of the non-additivities in PAWN, we computed second and third-order effects, shown in Fig. 4 (the second and third-order Sobol’ indices for each function and parameter are shown in Figs. S11–S14). In the $\max\notin\theta$ setting, the interactions that have a significant effect on the model output involve the initial sample size ($N$) with the number of conditioning intervals ($n$) or the stochasticity in the sampling of $Y_{U}$ ($\varepsilon$). Such second-order effects might contribute up to 15% of the PAWN index uncertainty. This also applies to the \sayOptimum setting, where the interaction between $N$ and $\varepsilon$ or $n$ and $\varepsilon$ has an even higher effect (up to 25% uncertainty). These three design parameters have significant third-order interactions in all three settings (Fig. 4b).

Regarding $S_{Ti}^{*}$ (Fig. 3b), both $N$ and $\theta$ have a non-negligible first-order effect on the index. The wide boxplots suggest that both the extent of their influence and the degree of non-additivities in defining the final $S_{Ti}^{*}$ index depends on the model under scrutiny (Figs. S8–S9).

4 Discussion and conclusions

Sensitivity analysis is an important tool to check the robustness of a model in the context of its validation. However, also the measurement of the sensitivity of the output variables to input parameters rests on modelling hypotheses, i.e. sensitivity analysis is a modelling process per se, based on the use of an algorithm – like the models being investigated. Here we tested the robustness of this aspect of the modelling process by assessing the dependency of the PAWN index to its design parameters, and matched the results against the Sobol’ total-order index ($S_{Ti}^{*}$). Our work thus sheds some light on the uncertainties concealed in the last stage of any sensitivity analysis: the implementation of the sensitivity algorithm. This practice could help ensuring the adequacy and the range of applicability of the tools which modellers deploy to improve the quality of their modelling exercises.

Two elements emerge from our work: the PAWN index is more sensitive to the design parameters than $S_{Ti}^{*}$, and this sensitivity has a complex pattern which makes the use of PAWN, or better the tuning of the PAWN design parameters, a delicate task. The chances of incurring in false positives (mistaking a non-influential parameter for influential, with a waste of computational resources) and/or false negatives (mistaking an influential parameter for non-influential, with a loss of useful information) are significant even when the design parameter space does not include the maximum as a possible summary statistic, but only measures of central tendency. Even in an ideal setting, where the uncertainty in $n$ and $\theta$ is highly reduced and $N$ ensures convergence, the PAWN index might be incapable to differentiate between non-influential model inputs and influential model inputs whose effect in the model output is fully through interactions. [31] also observed that the high level of noise in PAWN might produce biased rankings when the model inputs have different orders of magnitude in their contribution to the response.

PAWN especially underperformed when used in the Sobol’ G [23] ($k=8$) and the [24] function ($k=20$). This raises a red flag for analysts when using PAWN in high-dimensional models such as those commonly being employed in the Earth and Environmental Sciences domain, which might easily include tens of parameters [32]. In these contexts, the use of PAWN might indeed allow to significantly reduce the number of model runs required to conduct a sensitivity analysis, but at the expense of significant risks of obtaining a biased result due to its structural design.

The fact that the sensitivity of PAWN to its design parameters is complex, including important interactions up to the third order effect, implies that finding the perfect range of design parameters to use PAWN safely and efficiently is not easy. The significant non-additivity of PAWN only unfolds once the values of its main design parameters are moved simultaneously within reasonable uncertainty ranges; this is, when all the forking paths and divergences leading towards its computation are assessed at once. Instead, in the paper where the \saygeneric approach is presented, [2] analysed the influence of ($N,n,\theta$) on PAWN by combining different discrete point-estimates for $N$ and $n$ (see their Fig. 3), or by changing the value for either $N$ or $n$ while keeping the other design parameters fixed (see their Figs. 6 and 7). This approach is very similar to a one-at-a-time (OAT) sensitivity analysis, a method that may fail to detect interactions between model parameters due to its incomplete examination of the uncertainty space [33]. Our work shows that even by increasing $N$, ensuring a \sayhigh $N/n$ ratio or fixing $\theta$ at a central tendency measure, significant interactions between $Nn$, $N\varepsilon$ or between $Nn\varepsilon$ make PAWN considerably difficult to tame.

In our study we also observed the existence of potential interactions between the total number of model runs $N$ and the estimator $\theta$ for $S_{Ti}^{*}$. However, the variance in $S_{Ti}^{*}$ deriving from the uncertainty in its design parameters is comparatively much smaller, thus reducing the chances of obtaining a biased ranking or producing a wrong screening. It should be stated that our research design for $S_{Ti}^{*}$ included the use of different estimators out of our willingness to fully explore its potential uncertainty space. However, there is actually no uncertainty in $\theta$ when computing Sobol’ total-order indices: [27] showed that the best estimate of the degree to which a given model input interacts with the rest is obtained with the [20] estimator, which should be the default choice. This means that the uncertainty in the design parameters of $S_{Ti}^{*}$ actually narrows down to the total number of model runs $N$ only. It is thus to be expected that given two algorithms, the one depending on a higher number of design parameters will be the more delicate to use. A similar discussion led [25] to conclude that the screening method of Morris could be effectively substituted by $S_{Ti}^{*}$.

We would like to stress that there are many other uncertainties beyond those embedded in sensitivity indices that have the potential to condition the results of a given sensitivity analysis. These do not only involve quantifiable elements such as the selection of the distributions, sample size or the model output form [10, 11, 16], but also uncertainties in the entire knowledge and model-generating processes [19]. Algorithm uncertainties ingrained in the design parameters of sensitivity indices, however, are especially important for they condition the computational path towards the final sensitivity value. Although several alternative routes might be available, these multiple options are eliminated once the computation is triggered, yielding the illusion of a deterministic process. We believe that a good way to know whether a sensitivity index is robust is to check how volatile it is when all paths involved in its effective computation are walked at once.

The present findings do not suggest discarding PAWN as a sensitivity measure. Moment independent measures have a role to play in sensitivity analysis of output with long-tailed distributions. Additionally, they may find an ideal use in settings where the output of interest is itself in the form of a difference between two cumulative distributions. The present analysis aims at encouraging developers of sensitivity indices to fully explore the structural uncertainty of their algorithms in order to deliver transparent and robust sensitivity tools.

Acknowledgements

We thank Francesca Pianosi, Razi Sheikholeslami, Thorsten Wagener and two anonymous reviewers for their constructive comments on previous versions of this manuscript. This work has been funded by the European Commission (Marie Skłodowska-Curie Global Fellowship, grant number 792178 to A.P.).