Permutation Tests

Permutation Tests Jaime Tur´ on Rodr´ ıguez Supervisor: Ana Mar´ ıa Mu˜ noz Reyes Bachelor Thesis Faculty of Mathematics University of Seville July 2021 To Josefa, as I have longed for you to be alive to witness this. 3 Contents Abstract 7 Resumen 9 1 Introduction 11 1.1 The metaphysical duality of non-parametric statistics . . . . . . . 12 1.2 Prerequisites.............................. 13 2 Theory of Permutation Tests 17 2.1 Conditional Aspects of Permutation Tests . . . . . . . . . . . . . 19 2.2 What is a Permutation Test . . . . . . . . . . . . . . . . . . . . . 23 2.2.1 Constructing Permutation Tests . . . . . . . . . . . . . . . 23 2.2.2 Analyzing the p-Value . . . . . . . . . . . . . . . . . . . . 26 2.3 Statistics Permutationally Equivalent . . . . . . . . . . . . . . . . 28 2.4 The Test Statistic and the Critical Region . . . . . . . . . . . . . 31 3 More Properties of Permutation Tests 35 3.1 Unbiasedness ............................. 36 3.1.1 Characteristics of Conditional Unbiasedness . . . . . . . . 38 3.1.2 Strictly Uniform Conditional Unbiasedness . . . . . . . . . 41 3.2 Conditional Power Function . . . . . . . . . . . . . . . . . . . . . 43 3.2.1 Let’s make it unconditional! . . . . . . . . . . . . . . . . . 46 3.3 Constructing Confidence Intervals . . . . . . . . . . . . . . . . . . 47 3.4 Asymptotic Properties . . . . . . . . . . . . . . . . . . . . . . . . 49 5 6 Contents 3.4.1 Consistency .......................... 49 3.4.2 Asymptotic Behaviour of the Critical Value . . . . . . . . 49 4 Applying Theory 53 4.1 Testing Fixed Effects . . . . . . . . . . . . . . . . . . . . . . . . . 53 4.1.1 Balancing data . . . . . . . . . . . . . . . . . . . . . . . . 58 4.2 TestingSymmetry .......................... 59 4.2.1 Asymptotic behaviour for symmetry statistics . . . . . . . 64 4.3 A study for Ordered Categorical Variables . . . . . . . . . . . . . 64 Bibliography 71 Abstract In this Bachelor dissertation permutation tests are reviewed both theoretically and practically. Chapter 1 is a brief introduction to set the reader on the starting point of this project. Non-parametric methods and some basic concepts would be defined. Chapter 2 is the theoretical core of this dissertation. Permutation techniques are exhaustively introduced, as well as the construction of tests based on this idea, the behaviour of p-values, a succinct review of the critical region and an overview of permutational equivalence. Chapter 3 complements the previous chapter as theoretical framework development continues. Test unbiasedness, conditional power function, confidence intervals and general asymptotic properties are illustrated here. Finally, Chapter 4 shows some practical experiments where theoretical concepts approached in Chapters 2 and 3 are applied to real-world scenarios. 7 Resumen En este trabajo se presentan los tests de permutaciones de forma tanto te´orica como pr´actica. El Cap´ıtulo 1 es una breve introducci´on para situar al lector en el punto de partida del proyecto. Se definen los m´etodos no param´etricos y algunos conceptos claves. El Cap´ıtulo 2 es el n´ucleo te´orico del ensayo. Las t´ecnicas de permutaciones se explican en profundidad, as´ı como la construcci´on de tests basados en esta idea, el comportamiento de su p-valor, una breve introducci´on de la regi´on cr´ıtica y el concepto de equivalencia permutacional. El Cap´ıtulo 3 complementa el cap´ıtulo anterior continuando el desarrollo te´orico. La insesgadez de los tests, la funcion potencia condicionada, intervalos de confianza y algunas propiedades asint´oticas generales son los protagonistas de esta parte. Finalmente, en el Cap´ıtulo 4 se exponen algunos experimentos pr´acticos donde los conceptos te´oricos analizados previamente en los cap´ıtulos 2 y 3 se aplican a casos reales. 9 Chapter 2 Theory of Permutation Tests In this chapter we will introduce the basic concepts of permutation tests. We will use the two-sample design as a guide, following the spirit of [7]. Conditionality and exchangeability would be our cornerstone to develop the theoretical body as described in 1.2. As it was also said in advance, the conditioning will be done with respect to a set of sufficient statistics. Let Xbe a random variable taking values on the sample space X, and P∈ P the underlying distribution belonging to a non-parametric family of distributions (if Pwere a parametric family, classical parametric methods could be employed). The existence of the density function of P,fP(x), is assumed. This function would sometimes be interpreted as the density function itself whereas it could be also used to refer to the likelihood function. No distinction will be made as the context will be self-explanatory to determine which one we are referring to. Let X1={X11, ..., X1n1}∈Xn1and X2={X21, ..., X2n2}∈Xn2be two independent identically distributed samples data from the models P1and P2respec- tively, both P1, P2∈ P. We will write X={X11, ..., X1n1, X21, ..., X2n2}∈Xn, uniting the two i.d.d sample data; whose model is P=Pn1 1·Pn2 2, as X11, ..., X1n1 are i.i.d from P1,X21, ..., X2n2are i.d.d. from P2, and the models P1and P2 17 18 present independence; and where n=n1+n2. From now on we will write X= (X1,X2) = {X(i), i = 1, ..., n;n1, n2}to lighten notation, where the first n1data come from the first sample, and Xn1+1, ..., Xnbelong to the second one. Finally, we will note X∗={X(u∗ i), i = 1, ..., n;n1, n2}the resulting data sample from applying the permutation u∗= (u∗ 1, ..., u∗ n) to X. Remark that when permuting mixing observations between both populations is permitted. In fact, this would be the key point of permutation tests. If P1 and P2were identically distributed, this exchange should not be of relevance. However, if it is not the case, that is, P16=P2, changing data would be relevant when measuring statistics on the permuted data set. We develop the theoretical body in this project standing on testing problems for stochastic dominance generated by some possible fixed effects as a result of the application of a treatment. Stochastic dominance implies one sided alternatives. Two sided alternatives would be discussed later. Under this assumption, our null hypothesis is H0:{P1=P2} whereas the alternative hypothesis is H1:{X1+δ1> X2+δ2} with δ1and δ2the respective fixed effects in each population. In case we had quantitative variables, the null hypothesis could equivalently be stated as H0:{F1(t) = F2(t)∀t∈ R} and the alternative as H1:{F1(t)≤F2(t)∀t∈ R} where the strict inequality is given in a set of non null probability for both distributions. It is important to underline that under the null hypothesis data are Chapter 2. Theory of Permutation Tests 19 exchangeable. Although it is an assumption, it naturally stems from the fact that we should not be able to distinguish whether data come from the first or from the second sample as both populations behave equally; and that under the alternative, due to stochastic dominance, it is obvious that CDFs do not cross. Without loss of generality, we can assume that δ1=δ > 0 and P[δ2= 0] = 1. This might be interpreted as the first sample receiving a real treatment, where we want to test whether the treatment leads or not to some effects, while the second sample is under a placebo effect. Now, our hypotheses can be written as H0:{δ= 0}vs. H1:{δ > 0} This alternative implies that δi>0 for at least one individual of the first sample. We will use X(δ) = {X11 +δ1, ..., X1n1+δn1, X21, ..., X2n2}to note the data set under the alternative, whereas X(0) denotes the data set in H0(with no effects produced). It is also important to remark that the alternative hypothesis H1:{δ > 0} does not imply δi>0∀i= 1, .., n1, so the treatment may have some effects on some subjects of the first simple whereas it might not produce any effects on others. 2.1 Conditional Aspects of Permutation Tests We now introduce the first concept related to permutation tests: the conditional reference space. Definition 2.1.1. The conditional reference space is the set of points of the sample space Xwhich carry the same information as Xin terms of the likelihood function. It is noted X|X According to the definition, this space contains the points X∗such that the ratio fP(X)/fP(X∗) does not depend on P. Under the null hypothesis the density function fP(X) = Qj,i P(Xji) is assumed to be exchangeable in its arguments as 20 2.1. Conditional Aspects of Permutation Tests fP(X) = fP(X∗) and thus all permutations X∗of Xare in X|X. That is, the conditional reference space can be expressed as X|X={Su∗X(u∗ i), i = 1, ..., n}, the set of all possible permutations of the observed data set X. One important consequence is that every element X∗∈ X|Xis a set of sufficient statistics for the distribution Pin H0. Another consequence is that the sample space Xis partitioned into orbits X|X so that any point X∈ X belongs to one and only one orbit. We then have two disjunctive cases: if X1∈ X|X2then the orbits are equal, X|X1=X|X2, whereas X2/∈ X|X1instantly implies that orbits of X1and X2do not intersect, that is, X|X1TX|X2=∅. These consequences allow us to think about permutations as one element when working in the permutation approach, which is naturally consistent. We would not difference between points in the same orbit, but between orbits. Attending to this and enhancing the idea, the sample space could be written as X=SX/X|X where X/X|Xdenotes the partition of the sample space into orbits (the sample space is perceived as a quotient space, which embodies the notion of conditional reference spaces acting as units). Conditional reference spaces X|Xare also known as permutation sample spaces. The next aspect related to conditional procedures is conditioning with respect to a set of sufficient statistics in H0. When we condition with respect to X|Xthe null conditional probability of an event A∈ A given X|Xis independent of Pdue to its sufficiency. We then have that the permutation distribution induced by any statistic T:X → R, noted FT(t|X|X) = FT∗(t) = P[T∗≤t|X|X] does not depend of P, so it is P-invariant. This is a great advantage as any conditional inference is free of any distribution belonging to nonparametric families. Furthermore, we can define and give a closed form of the permutation probability of every A∈ A for finite sample sizes thanks to the finiteness of points in X|X, whose cardinal number is M=PX|X I(X∗∈ X|X)<∞, and to the independence on P: Chapter 2. Theory of Permutation Tests 21 P[X∗∈A|X|X] = PX∗∈AfP(X∗) PX∗∈X|XfP(X∗)=PX|X I(X∗∈A) M because for every permutation X∗∈ X|Xwe have fP(X) = fP(X∗). We can thus define the permutation measurable space (X|X,A|X), where A|X=ATX|X(the projection of events on the permutation sample space). Attending to the permutation probability introduced before, the following proposition is stated. Proposition 2.1.1. Let us assume that PX|X I(X∗=x) = 1 if x∈ X|Xand 0 elsewhere. Then permutations X∗are equally likely in H0: P[X=x|X|X] =    1/M if x∈ X|X 0 if x/∈ X|X and the elements of the same orbit X|Xare conditionally uniformly distributed over it. This proposition allows us to conclude that under the null hypothesis, the permutation distribution P[X∗=x|X|X] just depends on the resulting data set X. The data set Xcan be interpreted as the n-dimensional parameter for the permutation cumulative distribution function FT∗. Here we have an important difference in the alternative. A set of sufficient statistics is (X1,X2) being the data exchangeable within but not between samples, and so the observed data Xis not uniformly distributed over X|Xconditionally. We now introduce the empirical probability measure (EPM) and its analogous in case of quantitative variables, the empirical distribution function (EDF). Definition 2.1.2. Let X∗∈ X|Xbe a permutation of the observed data set X. The empirical probability measure of an event A∈ A is defined as: ˆ PX∗(A) = X i≤n I(X∗ i∈A)/n =X i≤n I(Xi∈A)/n =ˆ PX(A) (2.1) 22 2.1. Conditional Aspects of Permutation Tests which is a permutation invariant function. In case we have quantitative variables, we can also definde the empirical distribution function (EDF) as: ˆ FX∗(t) = X i≤n I(X∗ i≤t)/n =X i≤n I(Xi≤t)/n =ˆ FX(t) (2.2) which is permutationally invariant as well . An interesting result related to these concepts is stated. Although the EPM is the general version of the empirical function, we will deal with quantitative variables, so we state it just for the EDF. Proposition 2.1.2. Let X∈ X. Then ˆ FXis a permutation invariant function which characterises X|X. Thus, the conditional sample space can be defined as the set of points in Xwhich share the same EDF. This proposition is consequence of the partition of Xconformed by the different orbits X|X. A consequence of this result is that the EDF is a sufficient statistical function for Pin the null hypothesis. Hence, conditioning on the conditional reference space is equivalent to doing so on ˆ FX. Another consequence is that for any statistic T, P[T(X∗)≤t|X|X] = P[T(X∗)≤t|ˆ FX] = P[T(X∗)≤t|ˆ FX∗] for all t∈Rand all permutation X∗∈ X|X. According to this equivalence, the permutation null distribution of Tcan be interpreted as a process of without replacement random experiment from a uniformly distributed population whose distribution is the EDF. This interpretation is clearly different to that of bootstrap methods, where the selection is made with replacement. An asymptotic result with mild assumptions is obtained for the null permutation distribution of T. If we have sufficiently large n1and n2then the null permutation distribution P[T(X)≤t|ˆ FX] approximates its unconditional counterpart P[T(X)≤t] Chapter 2. Theory of Permutation Tests 23 2.2 What is a Permutation Test Let us suppose a test statistic Tfor which, without loss of generality, large values of Tare evidence against H0. Of course, if it is not the case, we can consider a suitable transformation of the test to be in this situation. Definition 2.2.1. On the previous conditions, the permutation support induced by (T,X) is defined as the set TX={T(X∗) : X∗∈ X|X} This set is essentially conformed by all the possible values which the statistic Tassumes by permuting the data set X. We now assume the null hypothesis as true. In light of 2.1.1, X∗is uniformly distributed over the conditional space. Let us put Mas the cardinal number of X|X∗, and let us consider T∗ (1) ≤T∗ (2) ≤ · · · ≤ T∗ (M), the members of TXin a non-decreasing order. If we now fix a value α∈(0,1), Tα(X) = Tα=T∗ (Mα) denotes the permutation critical value related to the statistic Tand the data set X, where Mα=PX|X I[T(X∗)< Tα] is the number of values assumed by Tin the permutation support lower than Tα. Note that the behaviour of the critical value not only depends on the data set observed Xbut also on X|X, as Mαwould vary depending on the values assumed by Tin the respective permutation sample space of X. That is, if we get another observed data set X∗verifying X∗∈ X|X, the critical value satisfies Tα=Tα(X) = Tα(X∗) as the two orbits coincide and so do the permutation supports. The critical value is orbit-invariant, so for each α∈(0,1), Tαis a fixed value in TXwhich varies as Xvaries in Xjust in case we jump into a different orbit. 2.2.1 Constructing Permutation Tests We now move on to construct the tests based on permutational aspects. As it is widely known, there are two kind of tests: randomized and non-randomized tests. Randomized Permutation Tests 24 2.2. What is a Permutation Test We first focus on randomized tests. Let us note φRthe test associated with the pair (T, X): φR=       1 if T0> Tα γif T0=Tα 0 if T0< Tα where T0=T(X) is the value of the test statistic Texamined on the observed data set Xand γis a probability given by γ=α−P[T0> Tα|X|X] P[T0=Tα|X|X] The test itself is not complete. When T0=Tα, we should define an auxiliary rule of decision based on a independent experiment. This can be made by running a variable U∼U(0,1) and rejecting H0in case U≥γ. Let us write the conditional expectation in H0: E[φR(X)|X|X] = P[T0> Tα|X|X] + γ·P[T0=Tα|X|X] =P[T0> Tα|X|X] + [α−P[T0> Tα|X|X]] ·P[T0=Tα|X|X] P[T0=Tα|X|X] =α what shows that randomized permutation tests are exact tests for all X∈ X and any α∈(0,1). Indeed, they verify a stronger condition, the uniform similarity property: Proposition 2.2.1. Under exchangeability of data X, the conditional rejection probability of a randomized test φRis X-P-invariant in H0for all P∈ P and all X∈ X. The proof can be found in [5]. It can be appreciated once again how robust permutation test are in a nonparametric framework, freeing not only the distribution but also the test and its critical region from any underlying unknown population and observed data. Chapter 2. Theory of Permutation Tests 25 Non-Randomized Permutation Tests Let us study now the other type of tests, non-randomized ones. These tests have the form: φ=   1 if T0≥Tα 0 if T0< Tα In contrast with the randomized version of permutation tests, the conditional expectation in H0(and coinciding with the type I error rate in non-randomized tests) is: E[φ(X)|XX] = P[T0≥Tα|X|X] = X X|X I[T(X∗)≥Tα]/M =αa≥α Considering the significance level function LX(t) = P[T∗≥t|X|X], and given a pair (T, X), the possible α-values for the test, called attainable α-values, are those where the significance level function is discontinuous, that is, when lX(t)>0, with lX(t) the derivative of LX(t). Thus, the set ΛX={LX(t) : lX(t)>0} (step points of the significance level function) contains these attainable α-values. This is always a discrete set which depends on T,Xand n. Because of this, non-randomized permutation tests do not admit all values of type I error rate in practice. Let us study the behaviour of ΛX. In case that it exists a constant c≥1 verifying that for all values in the permutation support PX|X I(T∗=t) = c(i.e., ΛXdoes not contain multiple points or these points have constant multiplicity), ΛXbecomes X-invariant, though it does not become free of its dependence on n. In this case, the set ΛXhas the form ΛX={mc/M, m = 1, ..., M/c}, and αavalues have constant jumps of c/M. Due to this, non-randomized tests become conservative when we choose a desired type I error rate αdwhich does not belong to ΛXand an attainable α-value with αd≥αa. Of course, if the desired rate is in the set of attainable α-values, exactness is obtained. 32 2.4. The Test Statistic and the Critical Region This lemma shows that permutation tests are better as any arbitrary test when testing invariance. The second result is related to how choosing the most powerful tests when testing invariance hypothesis again. Lemma 2.4.2. Let H0be an invariance hypothesis. Let fH1be the density of Pin H1, and X∈ X. We order the points of the sample space in terms of fH1 resulting in fH1(X∗ (1))≥fH1(X∗ (2))≥ · · · ≥ fH1(X∗ (M)). Then the most powerful test of size αis φR=       1 if fH1> fα γif fH1=fα 0 if fH1< fα where for any fixed α∈(0,1),fαis the critical value, Mαis the number of points which lay within the critical region and γ=α−P[f(X)> fα(X)] P[f(X)] =fα(X) With this result, along with the previous one, we can think that we do have solved the problem of choosing a suitable test. However, the disadvantage of this result is clear: we need to know the analytical form of the distribution P, and since we are working in a nonparametric framework, this results impracticable. We need to make a different approach to justify the selection of our test statistics. In fact, this selection is mainly heuristic. As Pis assumed to be unknown in the non-parametric framework, we can not make use of the previous results to establish a best statistic. In spite of this, there is an analogy in permutation tests and the parametric solution when we are working with large sample sizes. As stated in [6], Proposition 2.4.1. Let Tbe a best statistic for the parametric family P, and assume that the unconditional critical region does not depend on any specific Chapter 2. Theory of Permutation Tests 33 alternative. Then its permutation counterpart is asymptotically equivalent (and the permutation version of Tis asymptotically a best statistic for P). Though it is a powerful result as we are able to determine a best (asymptotically) permutation statistic for the family of distributions, this is a very particular case where Pis parametric. Since we have introduced permutation tests to work in a non-parametric framework, following the nature of the previous results to determine an optimal tests, this result is very restrictive. In light of these previous results, the choice of divergence turns out to be pretty heuristic. Results claiming optimal properties for tests are so rigid that the hypotheses are scarcely verified. What is more, permutation tests are presented to cover nonparametric situations, so meeting the requirements to apply some optimality results seems impractical. The other question to analyse is why large values are evidence against H0, a fact which we have arbitrarily imposed. As well as the selection of T, it is mainly heuristic, though we have some reasons backing up this fact. The first one is that this choice is supported on a suitable transformation which provides equivalence of tests. We can transform our test into another one in which the alternative distribution stochastically dominates the null distribution. What is more, we would prove in 3.1 that the null distribution of a statistic Tis dominated by its alternative distribution in one-sided critical regions. Thanks to this, power functions of tests would be monotonic increasingly with respect to the effects we are testing (more on this later in 3.2). Chapter 3 More Properties of Permutation Tests We continue on the same conditions previously stated, that is, we assume a one-dimensional variable Xtaking values on sample space Xwith probability distribution P, and that in H1the distribution of X1is shifted by a constant quantity δwith respect to that of X2. We are in stochastic dominance scenario, with F1(x)≤F2(x)∀x∈R. Moreover, we continue assuming that the constant effect δis positive. We would write the variable as X(ˆ δ) = {Xji =ˆ δji +Zji, i = 1, ..., nj, j = 1,2}in the alternative from now on, where Zji are random deviates exchangeable which follow an unknown distribution P. We also assume ˆ δ2i= 0 for i= 1, ..., n2 as we previously did, so there are no effects on the second group (recall the placebo effect). Under this assumptions, we can represent the data set as X(δ) = (Z1+δ, Z2), with δ= (δ11, ..., δ1n1). This notation emphasizes the exclusively activity of the treatment on the first sample, with no effects on the second one. It is important to remark that δ= (δ11, ..., δ1n1) is such that δ1i>0 for at least one i= 1, .., n1, and thus the null and the alternative state different situations. 35 36 3.1. Unbiasedness In these conditions, out hypotheses would be H0:{δ= 0}vs. H1:{δ > 0} 3.1 Unbiasedness Firstly let us introduce the concept of unbiasedness for a test statistic. We recover this definition from [5]. Definition 3.1.1. Let φbe a test. φis said to be unbiased if it verifies that βφ(θ0)≤α≤βφ(θ1) for every size α∈(0,1) and every specific alternative in H1, where βφ(θ)is the power function of the test. The notion of test unbiasedness is related to a comparison between the rejection behaviour of the test under the null and the alternative hypothesis respectively. It essentially states that the probability of rejecting is always higher when the alternative hypothesis is true than when the null hypothesis is true. Due to our determination of the rejection region for permutation tests (which is of the form [Tα,∞) as large values of Twere evidence against the null hypothesis), this definition can be equivalently stated as P[T(X)≥Tα|H0]≤α≤P[T(X(δ)) ≥Tα|H1] Owing to our test φis based on the test statistic T, we will say that Tpossesses unbiasedness instead of the test itself. We should not confuse this with the notion of unbiasedness related to the expected value of T. However, as the latter concept is defined in the parametric framework, there would be no chance of provoking confusion since we are working in a nonparametric scenario. We have a sufficient condition for determining the unbiasedness of T. Chapter 3. More Properties of Permutation Tests 37 Proposition 3.1.1. Let Tbe a test statistic. If the null distribution of T,TH0, is dominated by every distribution from H1, TH1, then the test φbased on Tis unbiased. Proof. The proof is straightforward from the definition of stochastic dominance. Let TH1dominate TH0for any alternative, and let us fix a type I error rate α∈(0,1). Then we have P[TH0≥Tα]≤α≤P[TH1≥Tα] However, this concept has been defined in an unconditional scenery. Since we are developing permutation tests from a conditional standpoint, we introduce the conditional counterpart of this notion. Definition 3.1.2. Let Tbe a test statistic. Tis said to be conditionally or permutationally unbiased if the p-values verify P[λT(X(δ)) ≤αa|X|X(δ)]≥P[λT(X(0)) ≤αa|X|X(0)] = αa for every X∈ X and any δ∈H1, where X(δ)=(Z1+δ, Z2),X(0) = (Z1,Z2), and αais any attainable α-value. Similarly to the classic definition of test unbiasedness, the interpretation of the definition is basically that rejecting is more likely when the alternative hypothesis is true than doing so when it is the null hypothesis the one which is true. In a similar way to the unconditional case, we can introduce a sufficient condition to obtain conditional unbiasedness. Proposition 3.1.2. Let Tbe a test statistic. If Tis such that λT(X(δ)) = P[T(X∗(δ)) ≥T0(δ)|X|X(δ)] ≤P[T(X∗(0)) ≥T0(0)|X|X(0)] = λT(X(0)) then Tis conditionally unbiased. 38 3.1. Unbiasedness We just have to use the monotony of the distribution function to get the result. Note that this sufficient condition just involves the respective p-values under both hypotheses. If we have smaller p-values in the alternative hypothesis, it would be more likely to reject using the same levels of significance than in the null. Of course, the reader could instinctively think that conditional unbiasedness is a more restrictive notion than the traditional one. This is in fact true, as we state in this proposition. Proposition 3.1.3. If Tis conditionally unbiased, then it is also unconditionally unbiased. This result naturally stems from the definition of conditional probability, as the unconditional rejection probability of the permutation test based on Tis ZX P[λT(X(δ)) ≤α|X|X(δ)]fP(X(δ)) ≥α because under the assumption that Tis conditionally unbiased, it is verified that P[λT(X(δ)) ≤α|X|X(δ)]≥αand RXfP= 1 as it is a density function. The converse, however, is not true. This reinforces the idea of conditional unbiasedness being a more stringent concept than the unconditional one. 3.1.1 Characteristics of Conditional Unbiasedness Let us dig into the concept previously introduced. Given an observed data set X, the two observed values of Tare T0(0) = T0(X(0)) (under H0) and T0(δ) =T0(X(δ)) (under H1). Consider a permutation of the observed data set, resulting into two different values of T,T∗(0) = T(X∗(0)) and T∗(δ) = T(X∗(δ)) respectively. We now define the increment of Tdue to the fixed effects δas Chapter 3. More Properties of Permutation Tests 39 ∆T(X(δ)) = T(X(δ)) −T(X(0)), and the difference of values resulted from permuting as ∆T(X∗(δ)) = T(X∗(δ)) −T(X∗(0)). These increments allow us to compare the behaviour of Tunder both hypotheses. Before continuing, we should make some assumptions on the statistic Twe are working with. The first one is that Twill have the form T=S1(X1)−S2(X2) (3.1) We assume Sjto be symmetric functions, which implies invariance under rearrangement of data (Sj(Xj) = Sj(X∗ j), with X∗ jany rearrangement of the data set Xj). It is typical to consider Sjsome kind of sample mean or median. The second assumption will be the non-decreasing monotony of Sj, so Sj(X+X0)≥Sj(X) for any data set Xand any non-negative X0≥0. This second assumption enhances the idea of large values of Tbeing evidence against H0as effects would be such that δ > 0. A final observation is that Tis non-increasing in its second n2arguments. Back on the increments, we have that ∆T(X(δ)) = T(X(δ)) −T(X(0)) ≥0, because Tis not decreasing in its first n1arguments along with the observation that it is non-increasing in its second n2arguments. We can also stablish a one-to-one pointwise relationship between the conditional permutation spaces of H0and H1. If we fix a permutation u∗= (u∗ 1, ..., u∗ n) of u= (u1, ..., un), then X∗(0) = {Z(u∗ i), i = 1, ..., n;n1, n2}gives a unique point in the alternative sample space, which is X∗(δ) = {Z(u∗ i) + δ(u∗ i), i = 1, ..., n;n1, n2}. It is also true the other way. Therefore, p-values in H0are λT(X(0)) = P[T(X∗(0)) ≥T0(0)|X|X(0)] whereas in H1, they can be expressed as 40 3.1. Unbiasedness λT(X(δ)) = P[T(X∗(δ)) ≥T0(δ)|X|X(δ)] =P[T(X∗(0)) + ∆T(X∗(δ)) −∆T(X(δ)) ≥T(X(0))|X|X(0)] Two key points here. The first has been continuously repeated and it is that the effects δare only active on the first sample, i.e., on the first n1observations of X. The second one is that, due to the previous assumptions, the difference in any data permutation ∆(X∗(δ)) tends to be smaller than ∆(X(δ)) because some effects are exchanged between the first and the second sample and Tis assumed to be non-decreasing on its first n1arguments. Let us clarify this. The data set originally is (Z1+δ, Z2) = (Z11 +δ11, ..., Z1n1+δ1n1, Z21, ..., Z2n2) in the alternative. When we consider a permutation u∗= (u∗ 1, ..., u∗ n), the new data set is (Z∗ 1+δ∗,Z∗ 2) = u∗(Z11 +δ11, ..., Z1n1+δ1n1, Z21, ..., Z2n2)=(u∗ 11(Z11)+ u∗ 11(δ11), ..., u∗ 1n1(Z1n1) + u∗ 1n1(δ1n1), u∗ 21(Z21), ..., u∗ 2n2(Z2n2)), where some values of the first sample, which are larger than those of the second sample in the alternative due to positive fixed effects, could have been permuted to the second sample. Therefore, ∆(X∗(δ)) −∆(X(δ)) is expected to assume non-positive values. Depending on how this occurs, different types of conditional unbiasedness are presented: 1. If ∆(X∗(δ)) −∆(X(δ)) ≤0 for all X∗∈ X|Xand for all data sets X∈ X, then λT(X(0)) ≥λT(X(δ)). This is denominated strictly uniform conditional unbiasedness. 2. If ∆(X∗(δ)) −∆(X(δ)) ≤0 is verified in terms of permutation distribution instead of pointwise, then the p-values also verify λT(X(0)) ≥λT(X(δ)) distributionally talking. 3. If for some data set X∈ X and some permutation X∗∈ X|X∆(X∗(δ)) − ∆(X(δ)) >0 is given, then we may not have conditional unbiasedness. This could happen when the two CDF cross. Chapter 3. More Properties of Permutation Tests 41 We will make use of the first kind of conditional unbiasedness, the strict uniform one, which is the most practical one. 3.1.2 Strictly Uniform Conditional Unbiasedness Let us prove strictly uniform conditional unbiasedness for test statistics of the form we have proposed earlier, T(X) = S1(X1)−S2(X2). The observed value in H0is T0(0) = S1(Z1)−S2(Z2) whereas, under H1, the statistic Tassumes the value T0(δ) = S1(Z1+δ1)−S2(Z2) Permutating the data set under H0and H1, we obtain T∗(0) = S1(Z∗ 1)−S2(Z∗ 2) and T∗(δ) = S1(Z∗ 1+δ∗ 1)−S2(Z∗ 2+δ∗ 2) = T∗(0) + ∆S(Z∗ 1, δ∗ 1)−∆S(Z∗ 2, δ∗ 2) We introduce the term δ2as some active effects of the first n1subjects may have been permuted to the second sample. Now let us make the following appreciations:  ∆S(Z∗ 2, δ∗ 2)≥∆S(Z∗ 2,0) = 0 = ∆S(Z2,0) as effects δ∗ 2are non-negative and Sjare invariant over rearrangement.  ∆S(Z∗ 1, δ∗ 1)≤∆S(Z∗ 1, δ1) as it may be some values of the effects produced on the second sample in δ∗ 1, which are null.  Finally, observe that ∆S(Z∗ 1, δ1) = ∆S(Z1, δ1) in distribution because the points of the conditional reference space are uniformly distributed. 48 3.3. Constructing Confidence Intervals we know this is an αsize test. If we subtract a certain quantity δ0for which the null hypothesis stated before is not rejected, then it is because that value clearly falls within the 1 −αconfidence interval for the real value of δ. Then based on this proposition, for a given statistic Tand an observed data set X, our aim is to determine the set of values for which the null hypothesis is not rejected. The fact of Tbeing permutationally equivalent to an estimator ˆ δ for δease this task (for example, T=PiX1i/n1−PiX2i/n2=¯ X1−¯ X2). Once again we resort to Monte Carlo simulation for finding the lower limit δlow (respectively, an algorithm for δup is presented later). 1. In first place, once ε(the desired width of the interval) and the level of confidence αare fixed, pick a negative number η < 0 and consider X1(η) = {X1i(δ)−(ˆ δ+η), i = 1, ..., n1} 2. Using also a Monte Carlo procedure based on Biterations, compute the empirical distribution function of T,ˆ F∗ B(T0 η), on T∗(η) = ¯ X∗ 1(η)−¯ X∗ 2(η) (with T0 η=T[X(η)]). 3. Repeat the previous steps varying the ηvalue until |1−ˆ F∗ B(T0 η)−α/2|< ε/2, resulting in δlow =ˆ δ+η. The algorithm for obtaining the upper limit is similar, just introducing a little variation. We will want to satisfy the condition |ˆ F∗ B(T0 η)−α/2|< ε/2, and δup =ˆ δ+ηwith η > 0 in this case. The initial selection of ηseems to be pretty heuristic. To make that approximation, we can attend to some values resulting from some parametric method if possible. In case the permutation support T(X) and the related attainable α-values are defined in discrete sets, the lower and upper limits are predetermined by the values in the permutation support whereas the probability of the interval is determined by the jumps of the permutation distribution function on the permutation Chapter 3. More Properties of Permutation Tests 49 support. These intervals will be of the form T(k)< t < T(k+1), leaving inside a probability of F([T(k+1)|X|X]) −F([T(k)|X|X]) 3.4 Asymptotic Properties 3.4.1 Consistency Keeping our aim of exploring different properties of permutation tests, we now move on to the notion of consistency, which is intimately related to the power function of a test. Recall a statistic Tis said to be consistent if, for fixed values of α∈(0,1) and δ > 0, the unconditional power function verifies lim n→∞ W(δ, α, n) = 1 This concept remarks that the test based on the test statistic Tis almost infallible when the sample size becomes extremely large, as we are providing an asymptotic probability of 1 for rejecting the null hypothesis when it is false. 3.4.2 Asymptotic Behaviour of the Critical Value In contrast with parametric tests, where the critical value is a fixed quantity, this now could vary as Xvaries in X, that is, they are random variables (recall 2.2.2). Let us study how does this variable Tαbehave when the sample size goes to infinity following [6]. Let us assume Xn∈ Xn,Mnand T(Xn), all of them considered successions. We shall also assume that Mn→ ∞ and PX|X(T∗< Tα)/Mn→1−αas n→ ∞ (that is, the test is of size α). Consider the randomized test φR=       1 if T(Xn)>(tα)n a(x) if T(Xn) = (tα)n 0 if T(Xn)<(tα)n 50 3.4. Asymptotic Properties with 0 ≤a(x)<1 and suppose that the following two conditions are satisfied: (1) There exists a constant lsuch that (Tα)n→lin probability. (2) There exists a function J(y) continuous at y=lsuch that for every yat which J(y) is continuous, P[T(Xn)≤y]→J(y) Under these assumptions, we have P[T(Xn)>(Tα)n|Xn]≤EX|X[φR(Xn)] ≤P[T(Xn)≥(Tα)n|Xn] which, thanks to (1) and (2), clearly leads to EX|X[φR(Xn)] →1−J(l) Both land J(y) depends on the sequence Pn. However, the dependence of J(y) is stronger than that of l. Due to this, we may introduce an alternative condition for (1): (1’) F(y|X|X)→F(y) in probability for every yat which F(y) is continuous, where F(y) is a distribution function, the equation F(y) = 1 −αhas a unique solution y=land F(y) is continuous at y=l. Actually this latter assumption implies (1). In effect, by definition P[(Tα)n≤y] = P[F(y|X|X)≥X X|X (T∗< Tα)/Mn] for every y∈R. Let ybe a point of continuity of F(y). As we have assumed, PX|X(T∗< Tα)/Mn→1−α=F(l), and y < l implies F(y)< F(l), so the right-hand side tends to 0 if y < l. It analogously tends to 1 in case y > l, and then (Tα)n→lin probability. What is more, let X∗ nbe a permutation of Xn, and X0 nanother element in the sample space, the three of them mutually independent. If T(X∗ n) and T(X0 n) have the limiting joint distribution F(y)·F(y0), then for every yat which F(y) Chapter 3. More Properties of Permutation Tests 51 is continuous F(y|X|X)→F(y) in probability and, if the equation F(y) = 1 −α has a unique solution y=l, it is also verified that (Tα)n→lin probability. We now define the test φ0: φ0=       1 if T(Xn)>(tα)n γnif T(Xn) = (tα)n 0 if T(Xn)<(tα)n assuming that the test has size α. If the previous assumption (1’) is satisfied, then we can firmly claim (tα)n→l. If (2) is also satisfied, we have E[φ0(Xn)] → 1−J(l). Considering the subfamily of distributions ˆ P ⊆ P for which (1’) and (2) holds, if ˆ Pcontains all sequences induced by H0, then φRand φ0are asymptotically equivalent in terms of power. An interesting application of this could be seen when testing symmetry. When approximating by a standard normal distribution the permutation distribution of a test statistic, (1’) holds. Because of this, the test based on the sample mean asymptotically behaves as the one-sided Student’s t-test of size α. Chapter 4 Applying Theory Once we have introduced the main theoretical aspects of permutation tests, we move on to illustrate them in some practical cases. The first case we will be studying is the one which embodies the theoretical framework that has been developed in this project. This is the hypothesis testing determining the existence of some kind of fixed effects on two populations versus the absence of these possible effects. The second case would present another hypothesis testing environment on which permutation approach would be useful. In this case we will be testing whether some distribution is symmetric or not. Finally, a third case would be introduced to show a possible permutation approach in a different nature problem, as we would be dealing with categorical variables. 4.1 Testing Fixed Effects As it was introduced, the first applied case concerning our study is testing the existence of non-null fixed effects. The dataset used for this task is related to the COVID-19 pandemic ([9]). This file contains 59 variables which gather information about the pandemic around the world, such as the date when data was collected, the country it belongs to; and other variables which are directly related to the disease itself by measuring some pandemic features such as daily new cases, total deaths, or the basic reproduction number of the disease (usually 53 54 4.1. Testing Fixed Effects noted R0). Our purpose here is to statistically determine, through a hypothesis testing, if vaccination is effective or not. As it is widely known, Israel was the fastest country when vaccinating its population. Because of these, we will isolate data coming from Israel and we will be working on this dataset. We select the daily new cases of this country from 02/21/2020 to 06/30/2021, studying the vaccination impact on this variable. Chapter 4. Applying Theory 55 Furthermore, we will make a second filter by splitting the data into two groups. The first one is the data within 02/21/2020 to 02/26/2021, while the second group is the data collected between 02/27/2021 and 06/30/2021. The reason behind this particular date is that it is the day when 50% of Israelis had already received at least one COVID-19 vaccine dose. We now study the behaviour of new daily cases in both groups. Let us denote X1= (X11, X12, ..., X1367) the observations coming from the first group and X2= (X21, X12, ..., X2128) the second group. Remark our sample size is n= 495, with n1= 367 and n2= 128. We have two populations, and we suspect the first one is affected by some kind of fixed effects due to the absence of vaccination. As there are less vaccinated inhabitants in the first group, the probability of being infected with COVID-19 could seem to be higher than the that in the second sample. This leads to our scheme (X1+δ, X2) and the hypothesis testing H0:{δ= 0}vs. H1:{δ > 0} First thing we should make is testing normality of both populations. Of course, if normality is not rejected, we could assume it and consequently apply traditional parametric methods such as an ANOVA or a t-student test. We run 56 4.1. Testing Fixed Effects the well-known Shapiro-Wilk test to this aim. As we see, both p-values are extremely small, leading to the rejection of normal behaviour in both populations. We should then look for some statistical alternatives, which would clearly be the permutation approach. The next step when testing is a pretty sensitive issue: we should determine what statistic Tshould we select to base our test on. As we have expounded previously, the difference of averaged means seems a suitable statistic to our aim. We thus select T=Pn1 i=1 X1i/n1−Pn2 i=1 X2i/n2. Observe that under the alternative Twould assume lower values when mixing observations between samples, as the first population is provided of not null side effects, whereas the second one Chapter 4. Applying Theory 57 assumes lower values. At this point, results explained in section 2.3 become extremely useful. We will work with the test statistic ˜ T=Pn1 i=1 X1i/n1due to its permutational equivalence to T. This fact eases the computational aspect of permutation tests. After setting the test up we should now compute the p-value. We go to the Monte Carlo algorithm presented in 2.2.2, fixing B= 1000. Figure 4.1: MC Algorithm Implementation following 2.2.2 The resulting p-value allows us to reject the null hypothesis and strongly conclude that δeffects are such that δ > 0. Running the algorithm for B= 10000, we obtain We can appreciate that the p-value decreases as the number of iterations B increases. Attending to the final comment made in 2.2.2, we can ensure through 64 4.3. A study for Ordered Categorical Variables Centring data with respect to the, we can appreciate data is slightly rightbiased, so higher values seem more probable than little ones. If the distribution were symmetric, we should expect a head as heavy as the tail. This is highly consistent with the decision we have made about the hypothesis testing. 4.2.1 Asymptotic behaviour for symmetry statistics Let us consider the test statistic T= n X i=1 Xi·Ui/( n X i=1 X2 i)1/2 The denominator is invariant at the permutation sample space. As E(T) = 0 and V ar(T) = 1, we can conclude Tis the sum of nindependent variables. Due to the PCLT, T∼ N(0,1) in the limit. 4.3 A study for Ordered Categorical Variables The third case study heavily differs from the two previous ones. Our aim now is to test different nature variables to those hitherto studied. We would be Chapter 4. Applying Theory 65 treating with categorical variables which would intrinsically possess some order notion. The current framework would be a categorical variable Xpartitioned into h≥2 classes, {Ak, k = 1, ..., h}. As we commented in advanced on the previous paragraph, the relationship Ak< Ajfor every pair where 1 ≤k < j ≤hwould be some self-explanatory issue. Given two independent random samples X1and X2our purpose here is to conclude whether categories are equally distributed in both samples or if distributions vary from one sample to another. In terms of hypotheses testing, this would be stated as H0:F1(k) = F2(k)∀k= 1, ..., h against the alternative hypothesis claiming some difference on the impact of at least one class, that is H1:F1(k)6=F2(k) for some k= 1, ..., h The notation Fj(k) = P{Xj≤Ak}, j = 1,2 is clearly consistent due to the order notion associated to classes. The dataset that we will be using in this study case is related to the number of deaths provoked by smoking in different countries periodically from 1990 to 2017 [10]. The socio-economic level of each country is also gathered in the dataset. This would be our target variable. We then want to study if the socio-economic level is a meaningful factor in smoking deaths. For this task, we select the number of passings in each country in 1992 and we compare it to that in 2017. To make a distinction, we would split data into two groups, the first group being countries where the number of deaths have decreased whereas the second groups is conformed by those where deaths have increased. This way we have two independent random samples to test the hypothesis testing of interest. 66 4.3. A study for Ordered Categorical Variables This split results in two samples: the first one is X1= (X11, ..., X156) whereas the second one is X2= (X21, ..., X2138). As it can be appreciated on the upper histogram, the majority of death differences accumulate around 0. Figure 4.3: A remarkable outlier could be appreciated. However, we can detect numerously outliers in the sample, including one which is specially far from the rest of observations. This case is China, a country where deaths has hugely gone up. Chapter 4. Applying Theory 67 Once we have split and skimmed the data, let us analyse the variable of interest. The socio-economic level of each country is specified in four different categories: low, lower-middle, upper-middle and high. As we previously explained, classes respect some kind of order relationship. Low Lower-middle Upper-middle High X1 7 17 8 24 X2 48 54 23 13 Table 4.1: Distribution of countries attending to socio-economic level As we can see, richer countries prevail among those where deaths have decreased (almost half of them are rich countries), whereas poorer countries tend to be in the second group, where deaths have gone up. As well as in the previous cases, we need now to determine a suitable statistic test to carry out the experiment. Attending to our aim of determining if categories are equally distributed over both samples, we propose a statistic involving the respective EDFs. This statistic is T= h−1 X k=1 [ˆ F1(k)−ˆ F2(k)]2·(ˆ F(h)[1 −ˆ F(h)])−1 where ˆ F1is the EDF of the first sample, ˆ F2the EDF of the second sample, and ˆ Fthe joint EDF. Observe that the denominator is permutationally invariant, so attending again to 2.3, we can dismiss this term and use the statistic T= h−1 X k=1 [ˆ F1(k)−ˆ F2(k)]2 which uniquely compares the squared differences of every class in each sample. 68 4.3. A study for Ordered Categorical Variables Remark that jumps are the probability of coming across with each class. In the first sample, the biggest gap is the one between the third and the fourth category, that is between upper-middle socio-economic level countries and rich countries. However, this changes in X2. The biggest jumps here are the one between 0 and the first category, poor countries, and the one between poor and lower-middle socio-economic level countries, meaning a higher presence of poorer countries in the second sample. Let us compute the p-value. We make use of the Monte Carlo algorithm (2.2.2) again to estimate it. Setting B= 1000 Chapter 4. Applying Theory 69 the resulting estimation is This p-value clearly leads to rejecting the null hypothesis. We can conclude that the socio-economic level of a country plays an meaningful role when it comes 70 4.3. A study for Ordered Categorical Variables to smoke habits, leading to more deaths in poorer countries. Bibliography [1] Fisher, R.A. (1935) The Design of Experiments 27-41 [2] Pitman, E.J.G. (1937) Significance Tests Which May be Applied to Samples From any Populations Wiley for the Royal Statistical Society. Supplement to the Journal of the Royal Statistical Society [3] Wasserman, L. (2006) All of nonparametric statistics Springer, Springer Texts in Statistics [4] Gibbons, J.D. and Chakraborti, S. (2003) Nonparametric statistical inference Marcel Dekker [5] Lehman, E.L. and Romano, J.P.(2008) Springer, Springer Texts in Statistics [6] Hoeffding, W. (1952) The large-sample power of tests based on permutations of observations The Annals of Mathematical Statistics [7] Pesarin, F. and Salmaso, L. (2010) Permutation Tests for Complex Data: Theory, Applications and Software John Wiley & Sons, Wiley Series in Probability and Statistics [8] Taleb, N.N. (2005) Fooled by randomness: the hidden role of chance in life and in the markets. New York: Random House 71 72 Bibliography [9] Ritchie, H. and more (2021). Israel: Coronavirus Pandemic Country Profile Our World In Data [10] PISA test score: Mean performance on the mathematics scale, 2000 to 2015 World Bank EdStats [11] Number of deaths from smoking in 1990 vs. 2017 Global Burden of Disease Collaborative Network