Assessing Simplifying Hypotheses in Density Estimation

PhD Dissertation Assessing simplifying hypotheses in density estimation Jose Ameijeiras Alonso ESCOLA DE DOUTORAMENTO EN CIENCIAS DOUTORADO EN ESTAT´ ISTICA E INVESTIGACI´ ON OPERATIVA SANTIAGO DE COMPOSTELA 2017 DECLARACI´ ON DO AUTOR DA TESE D. Jose Ameijeiras Alonso Alumno do Programa de Doutoramento en Estat´ıstica e Investigaci´on Operativa Como Autor da Tese de Doutoramento titulada “Assessing Simplifying Hypotheses in Density Estimation” Presento a mi˜na tese, seguindo o procedemento adecuado ao Regulamento, e declaro que: 1. A tese abarca os resultados da elaboraci´on do meu traballo. 2. No seu caso, na tese se fai referencia as colaboraci´ons que tivo este traballo. 3. A tese ´e a versi´on definitiva presentada para a s´ua defensa e coincide ca versi´on enviada en formato electr´onico. 4. Confirmo que a tese non incorre en ning´un tipo de plaxio de outros autores nin de traballos presentados por min para a obtenci´on de outros t´ıtulos. Asdo. Jose Ameijeiras Alonso AUTORIZACI´ ON DOS DIRECTORES DA TESE Dna. Rosa M. Crujeiras Casais Profesora do Departamento de Estat´ıstica, An´alise Matem´atica e Optimizaci´on D. Alberto Rodr´ıguez Casal Profesor do Departamento de Estat´ıstica, An´alise Matem´atica e Optimizaci´on Como Directores da Tese de Doutoramento titulada “Assessing Simplifying Hypotheses in Density Estimation” Presentada por D. Jose Ameijeiras Alonso Alumno do Programa de Doutoramento en Estat´ıstica e Investigaci´on Operativa Autorizan a presentaci´on da tese indicada, considerando que re´une os requisitos esixidos no artigo 34 do regulamento de Estudos de Doutoramento, e que como Directores da mesma non incurre nas causas de abstenci´on establecidas na lei 40/2015. Asdo. Rosa M. Crujeiras Casais Asdo. Alberto Rodr´ıguez Casal Agradecimientos Estos agradecimientos van por todas esas personas que me han ayudado a llegar hasta aqu´ı. En primer lugar, me gustar´ıa mostrar mi profunda gratitud a Rosa M. Crujeiras y a Alberto Rodr´ıguez Casal: sin su apoyo personal esta tesis nunca podr´ıa haber sido realizada. A ellos les debo la madurez investigadora que he adquirido a lo largo de estos a˜nos, la cual considero que es la parte m´as importante de haber realizado esta tesis. A Rosa le agradezco que siempre haya tenido un momento para m´ı y que me haya guiado de forma excelente durante estos a˜nos. A Alberto, el haber confiado en m´ı desde mis primeros pinitos en el mundo de la investigaci´on, el ejercer de psic´ologo cuando ha hecho falta y todas su cr´ıticas constructivas. I would like also to thank the collaborators that appear along the thesis, especially Christophe Ley and Francesco Lagona for hosting me in Ghent and Rome and for their wonderful guidance during my research stays. Como Christophe me dijo una vez: lo bueno que tiene este trabajo es que investigas con personas que se acaban convirtiendo en tus amigos. Cuando se est´a realizando una tesis tambi´en es importante cuidar la parte emocional. En este sentido, tengo que dar las gracias a Helena, por haber estado siempre ah´ı, aunque estuviese lejos en distancia, durante todos estos a˜nos, siempre ha sabido c´omo hacerme una compa˜n´ıa cercana en lo personal. A mis compa˜neros de comida, departamento, caf´e, a los itmatis and to my colleagues in Ghent and Rome. En particular a Mar´ıa, que a pesar de estar siempre ocupada, nunca ha dudado en sacar un momento para ayudarme cuando lo he necesitado, a Ar´ıs, ´ Angel, Edu y Juan Carlos por la compa˜n´ıa y por aguantar todas mis quejas, en especial, en esta parte final de la tesis. A mis compa˜neros de promoci´on (Mar´ıa, ´ Oscar y Xabi) y al resto de mis amigos de Vigo y de Santiago (´ Alex, Borjita, Dani, Jandros, Keko y Rub´en) que me han ayudado a desconectar cuando era necesario. Finalmente y ya que en esta tesis se trata el tema de la circularidad, en este caso, me gustar´ıa cerrar el c´ırculo de los agradecimientos con otra Rosa, mi madre, a la que le debo tanto. Sin ella nunca habr´ıa llegado tan lejos en la vida, muchas gracias por todos los esfuerzos que has realizado. Espero que esta tesis tambi´en sirva para rendir homenaje a los buenos momentos que he pasado con mis abuelos, Lola y Perfecto. Santiago de Compostela, septiembre de 2017. Jose Ameijeiras Alonso. Determinando Hip´oteses Simplificadoras na Estimaci´on da Densidade Resumo A an´alise estat´ıstica cl´asica de variables aleatorias unidimensionais soe enfocarse asumindo fen´omenos distribu´ıdos simetricamente arredor dun ´unico punto. Sen embargo, estes modelos poden non resultar adecuados cando tratan de capturar estruturas subxacentes m´ais complexas, frecuentes na an´alise de datos reais. Para resolver este problema, nos ´ultimos anos, dese˜n´aronse varias distribuci´ons de probabilidade que tratan de capturar aquelas situaci´ons nas que os datos presentan asimetr´ıa e multimodalidade (cando os datos se concentran arredor de m´ais dun valor). Esta problem´atica p´odese atopar tam´en en estruturas m´ais complexas de datos continuos, como nos datos circulares, ´e dicir, mostras que se poden representar como puntos na circunferencia dun c´ırculo unitario. Non obstante, antes de aplicar estes modelos m´ais flexibles, pero tam´en m´ais complexos, ´e importante determinar cando estes valen a pena. Neste senso, esta tese persegue dous obxectivos. Primeiro, dende un punto de vista non param´etrico, contrastar o n´umero de modas para datos lineares e circulares. Para iso, real´ızase unha revisi´on das diferentes propostas dispo˜nibles na literatura estat´ıstica e proponse un novo m´etodo que presenta un comportamento superior en escenarios tanto lineares como circulares. O segundo obxectivo, no contexto dos datos circulares, ´e determinar se a distribuci´on subxacente da mostra ´e (reflexivamente) sim´etrica arredor dunha direcci´on central. Respecto deste obxectivo, pres´entase un novo test e dem´ostrase que ´e ´optimo para contrastar a simetr´ıa circular. O comportamento das diferentes propostas, presentadas nesta tese, tam´en se analiza a trav´es de estudos de simulaci´on e son ilustradas empregando diferentes aplicaci´ons a datos reais. Palabras chave contraste de hip´oteses, estat´ıstica circular, estimaci´on non param´etrica, multimodalidade, simetr´ıa reflexiva Chapter 1 Introduction The classical statistical analysis of univariate continuous data is usually focused on the study of phenomena that are symmetrically distributed and concentrated around a single point. However, such models may fail to capture more complex underlying structures when applied to real data. To solve this issue, several distributions were designed trying to catch the situations where the data change when it is reflected around a central value (asymmetry) or when it is concentrated around more than one value (multimodality). On more complex structures of continuous data, such as circular data, that is, samples that can be represented as points on the circumference of a unit circle, this problem can be also found. In fact, one of the first references for statisticians in this field (Mardia, 1972, Ch. 1, compiling the first real progress in circular data analysis made during the 1950s and 1960s) deals with these topics from the very beginning, as asymmetric and multimodal distributions occur frequently in practice. However, before applying these more complicated models, it is important to determine whether it is worth it. With this objective, this essay is focused on assessing simplifying hypothesis in density estimation for linear and circular data. The goal of this thesis is twofold. First, testing the number of modes for linear and circular data and second determining if the underlying distribution of the data is (reflective) symmetric around a central direction in circular data. Section 1.1 is devoted to present the concept of mode and a brief overview on the nonparametric techniques for determining the number of modes is also given. Section 1.2 provides an introduction to the topic of circular data. The need of adapting the techniques for assessing the number of modes and some of the most widely used parametric distributions are revised in this section. In this section is also presented the circular (reflective) symmetry. A new testing procedure will be obtained using the Le Cam theory. In Section 1.3 some real datasets that motivate and illustrate the different 1 2CHAPTER 1. INTRODUCTION topics of this essay are described. Finally, in Section 1.4 the organization of this manuscript is provided and in Section 1.5 the principal contributions of this thesis are described. Contents 1.1 Assessment of the number of modes . . . . . . . . . . . . 2 1.2 Background on circular data . . . . . . . . . . . . . . . . . 7 1.2.1 Some circular distribution models . . . . . . . . . . . . . . . 8 1.2.2 Modes of circular data . . . . . . . . . . . . . . . . . . . . . 10 1.2.3 Circular reflective symmetry . . . . . . . . . . . . . . . . . 11 1.3 Real datasets . . . . . . . . . . . . . . . . . . . . . . . . . . 15 1.3.1 The 1872 Hidalgo stamp issue of Mexico . . . . . . . . . . . 15 1.3.2 Seropositivity in malaria eradication . . . . . . . . . . . . . 16 1.3.3 Date of the wildfires detected among the world . . . . . . . 17 1.3.4 Cracks in cemented femoral components . . . . . . . . . . . 18 1.4 Manuscript distribution . . . . . . . . . . . . . . . . . . . 19 1.5 Contributions of the thesis . . . . . . . . . . . . . . . . . . 21 1.5.1 Contributions on multimodality tests for linear data . . . . 21 1.5.2 Contributions on multimodality tests for circular data . . . 22 1.5.3 Contributions on symmetry tests in the circular setting . . 23 1.6 Acknowledgements . . . . . . . . . . . . . . . . . . . . . . . 23 1.1 Assessment of the number of modes The first problem to tackle is the determination of the number of modes. The concept of mode, for a continuous random variable X, is defined as the value (or values) at which the probability density function reaches a local maximum. If there is just one mode, then the density (or the associated distribution) is unimodal. Otherwise, it will be said that the density is multimodal (bimodal, if it has two modes; trimodal, if it has three; etc.). It will be seen that the local minima of the probability functions, namely antimodes, also play its role when the number of modes is tested. The local maxima and minima of the probability density function will be denoted as xiin this manuscript. Then, if fis a class–one function (i.e. it has a continuous first derivative) with jmodes and (j−1) antimodes, x1, . . . , x2j−1will denote the location 1.1. ASSESSMENT OF THE NUMBER OF MODES 3 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) Figure 1.1: Examples of linear probability density functions. From left to right models (as described in Appendix A): M4, M7, M18 and M21. The dashed lines represent the mode and the dotted lines the antimode locations. of these points. In Figure 1.1, some examples of unimodal and multimodal densities and their modes and antimodes are presented. The models considered for illustration are described in Appendix A and correspond to M4 (unimodal and symmetric), M7 (unimodal and assymetric), M18 (bimodal and symmetric) and M21 (trimodal and asymmetric). The identification of the (unknown) number of modes is quite common in applied fields, such as, geology, neurology, economics, ecology or astronomy. Some examples include the study of the percentage of silica in chondrite meteors (Good and Gaskins, 1980), the analysis of the macaques neurons when performing an attention–demanding task (Mitchell et al., 2007), the distribution of household incomes of the United Kingdom (Marron and Schmitz, 1992), the study of the body–size in endangered fishes (Olden et al., 2007) or the analysis of the velocity at which galaxies are moving away from ours (Roeder, 1990). In all these examples, identifying the number (and location) of local maxima of the density function (i.e. modes) is important per se, or as a previous step for applying other procedures. A first approach for determining the number of modes can be done in a parametric way, for example, using a mixture of distributions (a revision on this topic can be found in, for example, McLachlan and Peel, 2000). But, as it has been mentioned, the underlying structure of the data can be really complex, leading this parametric adjustment to a misspecification of the model. For that reason, a nonparametric perspective will be taken in this manuscript. Determining the number of modes can be done by estimating nonparametrically the probability density function f, as, in this case, it is not necessary to know anything about the real distribution of the sample. For doing that, the kernel density estimation can be used (see, for example, Wand and Jones, 1995, Ch. 2). Given a random sample X= (X1, . . . , Xn) of a random variable Xwith density function f, the kernel density estimation of fin a point xis 4CHAPTER 1. INTRODUCTION defined as: ˆ fh(x) = 1 n n X i=1 Kh(x−Xi),(1.1.1) where Khis the rescaled kernel function, that is, a unimodal and symmetric density function, centred at 0 and with bandwidth parameter h. An example of this estimator (and the one that it is going to be used later on) can be obtained using the Gaussian kernel: ˆ fh(x) = 1 nh√2π n X i=1 exp −1 2x−Xi h2!.(1.1.2) The bandwidth parameter plays a fundamental role in this estimation, as different number of modes can be obtained depending on the value of h. If the number of modes in ˆ fhis equal to k, then the location of the estimated modes and antimodes will be denoted as bx1,...,bx2k−1. In Figure 1.2 (left) the kernel density estimation with Gaussian kernel (1.1.2) obtained from a sample of size n= 50 from model M3 (see Appendix A) is represented, with bandwidths obtained from automatic rules, such us the rule of thumb and the Seather and Jones plug–in rule (see Wand and Jones, 1995, Ch.3). In this case, using the data–driven bandwidths, it is observed that the rule of thumb gives one mode and the plug–in rule shows four, while the real number of modes of fis j= 1. Although, in general, better results are achieved with the Seather and Jones plug–in bandwidth (for a discussion on this topic, see, for example, Wand and Jones, 1995, Ch.3), in this case it is the rule of thumb the one that detects correctly the underlying number of modes. The problem with these automatic bandwidths is that they are conceived to optimize the complete curve estimation (minimizing a global error criterion) and the problem here is focused only on the local maxima of f. Since the number of modes in ˆ fhis a monotone decreasing function of h, when the Gaussian kernel is used (see Silverman, 1981), a simple exploratory solution, for determining the number of modes, is representing this density estimation for different values of h as it is done in Figure 1.2 (right). In fact, this is the idea of some graphical tools, such as the mode tree and forest (see Minnotte and Scott, 1993; Minnotte et al., 1998) or the SiZer (see Chaudhuri and Marron, 1999). The problems with these exploratory tools are that an “expert eye” is needed to determine the number of modes and they do not provide a formal way for determining the number of modes (in the sense that no “statistical significance” is attached to the output result). Moreover, they do not enable a systematic study when the aim is to determine the number of modes in a large number of samples. For these reasons, a testing perspective will be taken. The lack of a test on the number of modes providing a satisfactory performance in practice motivated the development of Chapters 2 and 3, for linear and circular data, respectively. 1.1. ASSESSMENT OF THE NUMBER OF MODES 5 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x fh(x) −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x fh(x) Figure 1.2: Kernel density estimation with Gaussian kernel (1.1.2) obtained from a sample of size n= 50 (from model M3) represented with ticks in the bottom of the plot. Left (black lines): h= 0.116 (continuous; rule of thumb) and h= 0.0762 (dot–dashed; plug–in). Right (continuous lines from dark to light grey): critical bandwidths for one (h1= 0.110), two (h2= 0.096), three (h3= 0.087) and four (h4= 0.052) modes. Grey discontinuous lines: (approximately) estimated location of the modes (dashed) and antimodes (dotted). From a formal point of view, given a random variable with jmodes, the testing problem was formulated in the statistical literature (see, for example, Silverman, 1981) as: H0:j≤kvs. Ha:j > k, (1.1.3) where kis the number of modes to be tested. It should be noted that the null hypothesis, in the new proposals presented in this essay, is restricted to H0:j=k, but this change will not present relevant practical consequences as it will be discussed. There have been quite a few proposals in the statistical literature for solving (1.1.3) in the linear case and the different techniques can be classified in two groups: a first group of tests based on or using a critical bandwidth, and a second group of tests based on the excess mass. A brief description of the underlying ideas of these two concepts are given below, and more details will be provided in Chapter 2. The concept of critical bandwidth was introduced by Silverman (1981), using the aforementioned property of the number of modes being a monotone decreasing function of the density estimation (1.1.2). The critical bandwidth for kmodes, denoted as 6CHAPTER 1. INTRODUCTION hkis defined as the smallest bandwidth parameter for which ˆ fhhas kmodes. In the example shown in Figure 1.2 (right), the kernel density estimation with bandwidth h1(the darkest line) is the last estimation having one mode (near 0.5), for smaller bandwidth parameters a second mode emerges (near location x=−0.3). In a similar way, ˆ fh2has these two estimated modes and for smaller bandwidth parameters a third mode appears (near x= 0.9). Silverman (1981) uses this parameter as the test statistic and the null hypothesis (1.1.3) will be rejected when hkis too “large”, since this means that it is necessary to oversmooth the kernel density estimation to obtain ak–modal structure. A calibration of this test statistic when the support (where the modes are located) is bounded and known was provided by Hall and York (2001) in the case of testing unimodality (H0:j= 1). The critical bandwidth is also used by Fisher and Marron (2001) to construct a Cram´er–von Mises test statistic. In their paper, Fisher and Marron (2001) also provide a method for testing the number of modes when the random variable is circular using a U2of Watson (1961) test statistic. For the circular case, more details will be given in Chapter 3. The second group of procedures use as a test statistic the excess mass introduced by M¨uller and Sawitzki (1991). For testing kmodes, the underlying idea of this statistic is measure the difference in probability mass between considering that the density has (k+ 1) or kmodes. “Large” values of this test statistic will indicate that at least one more mode is present, i.e., that the null hypothesis is false. Although the first proposal in this group of tests, the one introduced by Hartigan and Hartigan (1985), does not use this test statistic, it will be presented in this block since the value of their test statistic (known as dip) is exactly a half of the value of the excess mass for the unimodal case. The dip statistic measures the maximum difference between the empirical distribution function and the integral of a “unimodal distribution function” (first convex, then with constant slope and in the final part concave) that minimizes that maximum difference. Both Hartigan and Hartigan (1985) and M¨uller and Sawitzki (1991) proposed the same way of calibrating the test statistic and although the dip and the excess mass are equivalent in the unimodal case, just the excess mass can be generalized for testing if the density is k–modal (with k > 1). The calibration of this test statistic for testing unimodality was provided by Cheng and Hall (1998). Finally, a comparative between the different methods is provided in Chapter 2. It will be seen that none of the existing proposals provides a satisfactory performance in practice. For that reason, a new procedure which combines the previous approaches (critical bandwidth and excess mass) is presented in Chapter 2 and compared with the existing methods, showing a superior behaviour (specially when testing H0:j=k, with k > 1). 1.2. BACKGROUND ON CIRCULAR DATA 7 1.2 Background on circular data As it was already mentioned, the second focus of this essay will be set on the circular case. The need of circular statistics appears for one dimensional data when the periodicity must be taken into account and the sample can be represented on the circumference. Examples of circular data appears when the goal is to model orientations or a periodic phenomenon with known period. These examples are frequent in many applied fields, some of them related with orientations: the orientation of the red wood ants in reaction to different stimuli (Jander, 1957); the flight orientation in pigeons returning home (Schmidt-Koenig, 1963); the wind direction at Gorleston (England) during 1968 (Fisher, 1995, Ch. 5); or waves directions in the Adriatic sea (Jona-Lasinio et al., 2012; Lagona et al., 2015). Other examples are related with periodic phenomena: the number of soldiers deaths (along each year) during the Crimean War (originally introduced by Nightingale in 1858 and revisited in Brasseur, 2005); or the frosting and defrosting times (during the day) in glaciers Monte Alvear and Vinciguerra in Ushuaia, Argentina (Oliveira et al., 2013). Having in mind that the data is periodic, in order to provide their characterization, two choices must be done: an initial direction and sense of rotation; but the results cannot depend on these choices. It is for that reason that, from an inferential perspective, circular data require a different treatment to that given to linear data. An example of this issue can be found in the sample mean. For example, if one has two directions in radians −π/6 and π/6, using as initial direction −π, then the sample mean is equal to zero; meanwhile if the initial direction is zero, then −π/6 becomes to 11π/6 and the sample mean is equal to π. For that reason, the sample mean direction must be defined taking into account the circular structure, given the circular sample in angles Θ= (Θ1, ..., Θn) from a circular random variable Θ, this mean can be defined as ¯ θ= arg (1 n n X j=1 cos Θj+i1 n n X j=1 sin Θj)(1.2.4) where iis the imaginary unit and arg denotes the function which returns the argument of a complex number. As it was already mentioned, a complete introduction on this topic can be found in Mardia (1972). More modern references include Fisher (1995), Mardia and Jupp (2000), Jammalamadaka and Sengupta (2001) or Ley and Verdebout (2017). A brief introduction to circular data is made in this essay to present notation and some special features of this kind of data. In this manuscript, the continuous and circular random variable Θ will be measured in radians an its support will be the interval [0,2π). To facilitate the explanation, an abuse of notation will be made denoting as falso to the circular density function, from now on and, also, in Chapters 3 8CHAPTER 1. INTRODUCTION and 4. Taking into account the aforementioned periodicity, the density function must satisfy the following conditions (Jammalamadaka and Sengupta, 2001, Ch. 2): •f(θ)≥0,∀θ∈[0,2π); •R2π 0f(θ)dθ = 1; •f(θ) = f(θ+ 2mπ),∀θ∈[0,2π) and ∀m∈Z. 1.2.1 Some circular distribution models Apart from the trivial case of the simple uniform distribution on [0,2π), some of the most well–known circular densities include the cardioid, the von Mises and the wrapped distributions such as the wrapped normal and the wrapped Cauchy (see, for example, Jammalamadaka and Sengupta, 2001). These four models (represented in Figure 1.3) are characterized by the concentration parameter, where “large” values of this parameter indicates when the distribution is more concentrated around zero. For simplicity, all the densities will be centred at 0, but a location parameter can be introduced to change the central direction, leading to densities of the form f(θ−µ), with central direction µ∈[0,2π). Cardioid distribution. The cardioid distribution centred at 0 with parameter of concentration ρ, namely Cρ, is a perturbation of the uniform circular density by a cosine function which has probability density function fCρ(θ) = (1/2π)(1 + 2ρcos θ),with θ∈[0,2π) and 0 < ρ < 1/2. Although in most of the references, introducing the circular topic, the concentration parameter is restricted to |ρ|<1/2, both here and in the following distributions, the support of the concentration parameter is taken to guarantee that this distribution is going to be unimodal. Being the mode at 0 or µ, when the cardioid distribution with location parameter µis considered, namely as C(µ, ρ) and with density fCρ(θ−µ). Von Mises distribution. The von Mises distribution, VM(µ, κ), was created to guarantee that the maximum likelihood estimator of µis equal to the circular sample mean defined in (1.2.4). The von Mises centred at 0, VMκ, has probability density function fVMκ(θ) = 1 2πI0(κ)eκcos(θ),with θ∈[0,2π) and κ > 0,(1.2.5) where the normalising constant contains I0which denotes the modified Bessel function of the first kind and order 0, I0= (1/2π)R2π 0exp(κcos(θ))dθ. 1.2. BACKGROUND ON CIRCULAR DATA 9 0.0 0.2 0.4 0.6 0.8 θ f(θ) 0π2π3π22π 0.0 0.2 0.4 0.6 0.8 θ f(θ) 0π2π3π22π 0.0 0.2 0.4 0.6 0.8 θ f(θ) 0π2π3π22π 0.0 0.2 0.4 0.6 0.8 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + 0 π 2 π 3π 2 + 0 π 2 π 3π 2 + 0 π 2 π 3π 2 + Figure 1.3: Examples of circular probability density functions. From left to right, linear (top) and circular (bottom) representation of the models centred on πwith different concentration parameters: von Mises, with κ= 0.1 (solid line), 1 (dashed line) and 5 (dotted line); cardioid, with ρ= 0.1 (solid line), 0.25 (dashed line) and 0.45 (dotted line); wrapped normal, with ρ= 0.2 (solid line), 0.5 (dashed line) and 0.9 (dotted line); and wrapped Cauchy, with ρ= 0.2 (solid line), 0.5 (dashed line) and 0.7 (dotted line). Wrapped distributions. Another way of obtaining circular distributions is wrapping around the circumference of unit radius a given linear distribution. That is, if Xis a random variable on the real line, the corresponding circular random variable Θ is given by Θ = Xmod 2π. Hence, if Xhas a linear density function fL, then the corresponding circular density function fCof Θ is fC(θ) = ∞ X p=−∞ fL(θ+ 2pπ). Some of the most well–known distributions of this kind are the wrapped Cauchy and the wrapped normal distributions. As its name indicates, the wrapped Cauchy distribution, fWCρ, is obtained wrapping the Cauchy distribution (see Johnson et al., 1995, Ch. 16), with location parameter 0 and scale parameter −log(ρ), where log is the natural logarithm. Its probability density function has the following expression fWCρ(θ) = 1 2π 1−ρ2 1 + ρ2−2ρcos θ,with θ∈[0,2π) and 0 <ρ<1. 16 CHAPTER 1. INTRODUCTION types and, in general, with a lack of quality control in manufactured paper, which led to important fluctuations in paper thickness, being thin stamps more likely to be produced than thick ones. Given that the price of any stamp depends on its scarcity, the thickness of the paper is crucial for determining its value. Although the watermark of some stamps (a mark indicating were the stamp was produced; see Figure 1.6, right) can help to catalogue a stamp issue, there is not a standard rule for classifying stamps according to their thickness (not being available such a classification in catalogues), becoming this problem even harder in stamp issues printed on a mixture of paper types with possible differences in their thickness. A stamp issue with this problematic is the Hidalgo catalogue (an stamp of this issue can be seen in Figure 1.6, left), from the nineteenth century, and this particular example has been shown in several references in the literature as a paradigm of the problem of determining the number of modes/groups. From a nonparametric point of view, some examples of its utilization for mode testing have been shown by Efron and Tibshirani (1994, Ch. 16), Izenman and Sommer (1988) or Fisher and Marron (2001). It has been also analysed using nonparametric exploratory tools by Wilson (1983), Minnotte and Scott (1993) and by Chaudhuri and Marron (1999). This dataset was obtained from Table 2 in Izenman and Sommer (1988). 1.3.2 Seropositivity in malaria eradication Malaria is a global health problem causing symptoms that typically include fever, fatigue, vomiting, headaches and, in the several cases, it can cause yellow skin, seizures, coma, or death. This disease, caused by Plasmodium parasites, is most commonly transmitted to humans by infected mosquitoes. Being two of the species of Plasmodium that can infect and be spread by humans, the Plasmodium falciparum (Pf) and the Plasmodium vivax (Pv). Due to the prevention measures in the last decades, in several countries are already targeting malaria elimination. A way of understanding the true malaria exposure in a given population is using the antibody data against different malaria antigens (see, for example, Corran et al., 2007, a test detecting the antigens is provided in Figure 1.7, left). The reason of using this data is that the antibody concentrations are directly correlated with the parasite exposure, thus, it provides information on current and recent infections. In seroepidemiological studies, one of the most popular antibodies is the Merozoite Surface Protein–1 (PfMSP1 when referring to Pf and PvMSP1 when it is associated to Pv). Employing the antigen data, when no more information about the serological status (seropositive or seronegative) of the individuals is provided, a way to determine if malaria is still in a given population is try to determine if there are more than one latent subpopulation. With the objective of seeing if malaria is in an elimination phase, data from Aneityum island 1.3. REAL DATASETS 17 Figure 1.6: Left: doce centavos stamp of the Hidalgo issue of Mexico printed on 1872 from Wikimedia Commons (2009b). Right: example of a watermark on the back of a Zululand (historical region of South Africa) stamp indicating that the stamp belongs to the Crown CA production, from Wikimedia Commons (2009a). in Vanuatu archipelago and Chabahar in southeastern Iran are used. In this two regions, where Pf and Pv parasities co–exist, the antibodies PfMSP1 and PvMSP1 were quantified. The complete description of these two dataset is given in Cook et al. (2010) and Zakeri et al. (2016). More complete information about the malaria topics and some references can be found in, for example, Sep´ulveda et al. (2015). These datasets were provided by Dr. Nuno Sep´ulveda from the London School of Hygiene and Tropical Medicine, United Kingdom, for a joint work with other coauthors on this topic (see Section 1.5.1). 1.3.3 Date of the wildfires detected among the world This dataset contains the location and date (in a daily resolution) of all the fires detected among the world by the MODerate resolution Imaging Spectroradiometer (MODIS), launched into Earth orbit by the National Aeronautics and Space Administration (NASA) on board the Terra (EOS AM ) and the Aqua (EOS PM ) satellites, from 10 July 2002 to 9 July 2012 (an example of some detected fires is included in Figure 1.8). Taking the date of all the fires detected in each cell of size 0.5◦(latitude and longitude), the goal is to assess if there is one season of fires along the year or more than one. This study can help to understand where humans are altering the fire seasonality causing more peaks of fires than those ones expected due to climatolo- 18 CHAPTER 1. INTRODUCTION Figure 1.7: Left: example of a dipstick detecting some specific antigens produced by malaria parasites in the blood of infected individuals, from Wikimedia Commons (2004a). Center: map of Iran from Wikimedia Commons (2004b), Chabahar is located in the South–East of Iran. Right: map with the Vanuatu archipelago from OpenStreetMap contributors (2017), marking (in orange and with a pointer) the Aneityum island. gical factors. In addition to the problem of testing the number of modes in the circle, another issue appears as the goal is to analyse the number of fire seasons globally. Then, due the territory division, the incorrect rejections of the null hypothesis caused the multiple testing must be controlled. In addition, the possible spatial correlation between the “closed” regions must be also taken into account. How to solve this problematic will be tackled in Chapter 3. This dataset was provided by Dr. Jos´e Miguel Cardoso Pereira and his research team from the Instituto Superior de Agronomia, Lisbon, Portugal, for a joint work with other coauthors on this topic (see Section 1.5.2). 1.3.4 Cracks in cemented femoral components This dataset was collected during an in vitro fatigue study of total hip replacements described in Mann et al. (2003) and also analysed in Oliveira (2014) for illustrating the CircSiZer maps (see Figure 1.9 center and right). The data from the study consist of the directions, measured in angles relative to the centre of the stem, of fatigue cracks around the cemented femoral components in hip implants (how the angles are obtained is shown in Figure 1.9). After a substantial stress cycle had been applied, each femur was sectioned in 10 mm intervals from the level of the implant collar to the distal tip of the stem. Measurements at 60 and 70 mm were not made because of features in the experimental setup. As a result, two groups of measurements were 1.4. MANUSCRIPT DISTRIBUTION 19 Figure 1.8: Fires detected (in red) by MODIS. Left: on board the Aqua satellite around Central African Republic on February 6, 2004, from NASA Visible Earth (2004). Right: on board the Terra satellite in the north–eastern part of Russia (near the Bering Sea) on July 18, 2003, from NASA Visible Earth (2003). obtained: those in the proximal (10–50 mm) region and those in the distal (80–110 mm) region. Mann et al. (2003) showed that the directions of the fatigue cracks are not uniformly distributed and that their distributions in the two regions differ. The objective in this case is investigate whether the cracks in the two regions are symmetrically distributed about some unknown centre. This dataset was provided by Dr. Kenneth A. Mann from the Upstate Medical University, New York, United States of America. 1.4 Manuscript distribution A brief summary of each chapter in this manuscript will be provided. The main topics will be highlighted, as well as the chapter distribution. In a subsequent section, the scientific contributions of this thesis will be presented. Chapter 2 contains a review of the different proposals for testing the number of modes. This task has been approached in the statistical literature from two different perspectives: the critical bandwidth and the excess mass. Since none of the existing proposals provides a satisfactory performance in practice, a new procedure 20 CHAPTER 1. INTRODUCTION ● Lateral Anterior Medial Posterior 10 20 30 40 50 ● Lateral Anterior Medial Posterior 10 20 30 Figure 1.9: Left: the data is obtained from the angle of the crack with respect to the lateral direction, taking counterclockwise as the positive sense of rotation, from Mann et al. (2003). Center and right: CircSiZer maps for kernel density estimates in proximal (central panel) and distal (right panel) regions; given a concentration parameter (as indicated in the radius) and an angle, blue colour indicates locations where the curve is significantly increasing, red colour shows where it is significantly decreasing and purple indicates where it is not significantly different from zero; thus, for a given concentration, a blue–red pattern indicate where there is a significance peak. which combines the previous approaches is presented and compared with the existing methods, showing a superior behaviour. The proposed technique is illustrated with the Hidalgo dataset and the example on the detection of seropositivity in malaria eradication. In Chapter 3 a new method for assessing the number of modes in the circular case is given and it is provided its behaviour in the finite sample case. The proposed method is used in a real dataset of wildfires in a global scale. The issues related with the multiple testing problem and the spatial correlation between observations in the different locations are also solved in this chapter. Chapter 4 is devoted to the development of a new procedure for testing circular reflective symmetry about an unknown central direction. The construction of the optimal (in the maximin sense) test will be done following the Le Cam theory. The asymptotic distribution of the test statistic under the null hypothesis of circular symmetry and under k–sine–skewed alternatives is provided. The behaviour of the test in the finite sample case is shown. Applications of the test in a dataset related with the cracks in cemented femoral components is also provided. Finally, the manuscript includes some final comments and discussion in Chapter 5 1.5. CONTRIBUTIONS OF THE THESIS 21 and four appendices. In Appendix A, the circular and linear models used in the simulation studies in Chapters 2, 3, 4 and with illustration purposes along the manuscript are defined and represented. Appendix B collects the technical proofs of Chapters 2 and 4. Appendix C includes some numerical approaches made in Chapter 2. In Appendix C, the utilities of a new Rpackage are presented. Package multimode includes different methods for testing and exploring the number of modes in a linear density function. In Appendix D, an introduction on Le Cam Theory is given. 1.5 Contributions of the thesis The main goal of this thesis is to propose new tests for assessing simplifying hypothesis in linear and circular density modelling. As mentioned in the Introduction, the most relevant contributions can be summed up in three parts: multimodality tests in the linear setting, multimodality tests for circular data and symmetry tests for circular densities. In the following paragraphs, a summary of the main contributions is presented. 1.5.1 Contributions on multimodality tests for linear data The objective is testing if a linear random variable has kmodes against more than kusing nonparametric methods, i.e, H0:j=kvs. H0:j > k, where jdenotes the unknown number of modes and k∈Z+. State of the art From a nonparametric perspective, the multimodality test were firstly introduced by Silverman (1981), Hartigan and Hartigan (1985), and M¨uller and Sawitzki (1991). Referring to the case when k= 1 (testing unimodality) just Cheng and Hall (1998) and Hall and York (2001, when the support of the random variable is known) provide asymptotically well calibrated procedures. In the more general case, when k > 1, just the proposals of Silverman (1981) and Fisher and Marron (2001) are available and none of them are well calibrated. A dataset in which analysing the number of modes is important is the known as the 1872 Hidalgo stamp issue of Mexico, introduced originally by Wilson (1983) and used as a motivational example since the analysis of Izenman and Sommer (1988), who employed the test proposed by Silverman (1981). Referring to the available software for testing the number of modes, up to the knowledge of the author, just the dip test of Hartigan and Hartigan (1985) and the critical bandwidth test of Silverman (1981) were implemented, see, for example, the diptest 22 CHAPTER 1. INTRODUCTION package of R(Maechler, 2015) or the silvtest function in Stata (Salgado-Ugarte et al., 1998). Contributions A new test outperforming, in the unimodal case, the previous asymptotically well calibrated procedures and allowing the extension to the more general case (when k > 1) is presented in Ameijeiras-Alonso et al. (2016). A new procedure which combines the previous approaches (smoothing and excess mass) is presented and compared with the existing methods. For illustration purposes, it is also included the 1872 Hidalgo stamp issue of Mexico example for analysing the number of modes with a well calibrated test. Another application of the proposed test is done in the malaria field for identifying seropositive groups in the population in Sep´ulveda et al. (2017). Finally, in Ameijeiras-Alonso et al. (2017a), the different methods, implemented in the multimode package of R, for testing and estimating the number and location of the modes are provided. 1.5.2 Contributions on multimodality tests for circular data The objective is testing the number of modes in a circular sample. Then, apply the method for determining, in each region of the world, if there is one season of fires or more lead to a multiple testing problem, where a spatially–adapted False Discovered Rate (FDR) correction method must be applied. State of the art For testing if the circular sample has kmodes nonparametrically, just the proposal of Fisher and Marron (2001) is available in the statistical literature and it seems to be not well calibrated. From a parametric point of view, determining the number of modes can be done using mixture models. In Benali et al. (2017), employing the mixture of two von Mises, the global map is catalogued in function of the number of modes. It will be relevant also for our example how to handle the multiple testing problem with spatially correlated data. Benjamini and Heller (2007) proposed an algorithm with this purpose, relying on prior information on the locations which may be the driving force of the dependency structure. Contributions In Ameijeiras-Alonso et al. (2017b) a new proposal for testing the number of modes 1.6. ACKNOWLEDGEMENTS 23 in the circular setting is provided showing a correct calibration in the finite sample cases. Its performance is compared with the Fisher and Marron (2001) proposal and then applied for testing if there is one season of fires in each region of the world. Then, taking into account the spatial dependence between the different locations, the results are corrected using an adaptation of the Benjamini and Heller (2007) proposal. 1.5.3 Contributions on symmetry tests in the circular setting The objective is testing circular reflective symmetry, i.e., if the density function of the data satisfies f(θ+µ) = f(−θ+µ), when the central direction µis unknown. State of the art For testing if the circular sample is reflective symmetric there are three proposals. The first one is the b2–bar test of Pewsey (2002), which allows to test the symmetry when the central direction is unknown, but there is no notion of optimality. When µ is previously given, both proposals of Pewsey (2004) and Ley and Verdebout (2014) can be used for testing circular symmetry, but just the second one is optimal (in the maximin sense) against the k–sine–skewed asymmetric models. Contributions In Ameijeiras-Alonso et al. (2017c) a new proposal for testing circular symmetry with unknown central location is given. The test is proved to be asymptotically well calibrated and optimal against k–sine–skewed alternatives. In addition, the asymptotic distributions of the test statistic under the null and, also, under these alternatives are given. The behaviour of this test for small sample sizes is analysed with an extensive simulation study. For illustration purposes, an example for studying if the cracks in the cemented femoral components of the hip implants are symmetrically distributed around a central direction is included. 1.6 Acknowledgements This work has been supported by the PhD grant BES–2014–071006, the Project MTM2013–41383–P (both from the Spanish Ministry of Economy, Industry and Competitiveness), the Project MTM2016–76969–P (Spanish State Research Agency, AEI, being both projects co–funded by the European Regional Development Fund, ERDF) and IAP network P7/06 StUDyS (Developing crucial Statistical methods for Under- 24 CHAPTER 1. INTRODUCTION standing major complex Dynamic Systems in natural, biomedical and social sciences) from the Belgian Science Policy. Part of the research done in Chapter 4 was carried out during a visit to Ghent University supported by the Grant EEBB–I–16–11503 (Spanish Ministry of Economy, Industry and Competitiveness). The support of the Mobility Grant EEBB–I-17–12716 (Spanish Ministry of Economy, Industry and Competitiveness) is also thanked. The Supercomputing Center of Galicia (CESGA) is acknowledged for providing the computational resources that allowed to run most of the simulations. Finally, Drs. Kenneth A. Mann, Jos´e M. C. Pereira and Nuno Sep´ulveda and their research teams are thanked for providing some of the datasets employed in this manuscript. Chapter 2 Mode testing for linear data Several proposals were developed in the statistical literature for testing the number of modes of a random variable in the linear setting. This chapter is devoted to review different proposals for testing if the underlying distribution of the data is unimodal or multimodal from a nonparametric point of view. Some of these tests will be adaptable to the general case, when the objective is testing if the underlying distribution is k– modal. As mentioned in Chapter 1, the existing proposals can be divided in two groups: a first group of tests based on or using a critical bandwidth, introduced by Silverman (1981), further studied by Hall and York (2001) and also used by Fisher and Marron (2001); and a second group of tests based on the excess mass, such as those ones proposed by Hartigan and Hartigan (1985), M¨uller and Sawitzki (1991) and Cheng and Hall (1998). These methods are compared in this chapter, showing the behaviour of these procedures under different scenarios. The main contribution of this chapter is a new testing procedure combining the use of a critical bandwidth and an excess mass statistic, outperforming the existing procedures, in testing unimodality and more general hypotheses. The illustration of the new technique for testing the number of modes will be done with two examples. The first one, in the philatelic field, introduced in Section 1.3.1, can be considered as a classical example for analysing the underlying number of modes. This real dataset, firstly studied by Wilson (1983), has become (since its use by Izenman and Sommer, 1988) a popular example, as determining the number of underlying groups, in this data, is not an easy task. This issue is shown in Figure 2.1, where observing the histograms or the kernel density estimations with different bandwidths (left panels) does not provide a clear answer about the number of modes. Not even the most well–known graphical tool for detecting significant features on a scale– space perspective, the SiZer map (Chaudhuri and Marron, 1999), provides a clear 25 32 CHAPTER 2. MODE TESTING FOR LINEAR DATA The methodology proposed by Cheng and Hall (1998) consists in generating samples from Ψ(·,ˆ b), where ˆ band the distribution family are chosen using ˆ d. The excess mass statistic given in (2.1.7) when k= 1, ∆∗ n,2, is computed from the resamples and, for a given significance level α, the null hypothesis is rejected if P(∆∗ n,2≤∆n,2|X)≥1−α. 2.1.3 A new proposal The previous tests show some limitations: first, just the proposals of Silverman (1981) and Fisher and Marron (2001) allow to test (1.1.3) for k > 1. Despite the efforts of Cheng and Hall (1998) and Hall and York (2001) for providing good calibration algorithms, it will be shown in Section 2.2 that the behaviour of all the proposals is far from satisfactory. Specifically, the test presented by Silverman (1981) is very conservative in general (although sometimes can show the opposite behaviour) and the proposal of Fisher and Marron (2001) does not have a good level accuracy. The new proposal tries to overcome these drawbacks by considering an excess mass statistic, as the one proposed by M¨uller and Sawitzki (1991) with bootstrap calibration. Unlike the method presented by Cheng and Hall (1998), a completely data-driven and nonparametric procedure will be designed, using the critical bandwidth under the null hypothesis to test H0:j=kvs. Ha:j > k, for a general k∈Z+. The proposal, in a nutshell. Consider then the testing problem (1.1.3) and take the excess mass statistic given in (2.1.7), under the null hypothesis. Given a sample X= (X1, ..., Xn), generate Bsamples X∗b(b= 1, . . . , B) of size nfrom (a modified version of) ˆ fhk. For a significance level α, the null hypothesis will be rejected if P(∆∗ n,k+1 ≤∆n,k+1|X)≥1−α, where ∆∗ n,k+1 is the excess mass statistic obtained from the generated samples. It should be also noted that the procedure can be easily adapted to handle Hall and York (2001) scenario: to test the null hypothesis that f has almost kmodes in the interior of a given closed interval I, if Iis known, use (a modified version of) ˆ fhHY,k to generate the samples. From this brief description, two questions arise: How is this modified version of ˆ fhkconstructed? Does the procedure guarantee a correct calibration of the test? In fact, the construction of the modified kernel density estimator, namely the calibration function and subsequently denoted by g, is intended to ensure the correct calibration, under some regularity conditions. Regularity conditions (RC1) The density function fhas continuous derivative. (RC2) There exist t1and t2, such that fis monotone in (−∞, t1) and in (t2,∞). (RC3) The unique points satisfying {x:f0(x) = 0 and f(x)6= 0}are the modes and 2.1. A REVIEW ON MULTIMODALITY TESTS 33 antimodes of f, denoted as xi, with i= 1,...,(2j−1); and f00(xi)6= 0. (RC4) f00 exists and is H¨older continuous in a neighbourhood of each xi. Although these premises are quite general, it should be noted that they exclude some classical distributions, such us, those ones having an unbounded density function (e.g., the Beta distributions with parameters equal to 0.5, see Johnson et al., 1995) or the uniform distribution. Under the regularity conditions, Cheng and Hall (1998) indicated that the distribution of the test statistic (2.1.7) is independent from the underlying density f, except for the values in the modes and antimodes: di=|f00(xi)| f3(xi)with i= 1,...,(2k−1).(2.1.8) Using the results of Cheng and Hall (1998), and assuming that fhas kmodes, the distribution of ∆n,k+1 can be approximated by ∆∗ n,k+1, just with bootstrap values obtained from samples coming from a calibration distribution with kmodes. This calibration distribution must satisfy that the estimated values b diconverge in probability to di, as n→ ∞, for i= 1,...,(2k−1). As aforementioned, it is for that reason that Cheng and Hall (1998) proposed using a parametric family depending on the value of d1when k= 1. Two issues appear related with their calibration procedure. First, it is not an easy task to construct a parametric family having the desired values of di, for i= 1,...,(2k−1), when k > 1. Second, as Cheng and Hall (1998) pointed out the second–order limiting properties of the test depend on the form of the density function. Then, it is expected a better behaviour if the calibration function is “more similar” to the real density function. Our method deals with these two issues to get a test having a good behaviour in the finite–sample case and allowing solving the general problem of testing kmodes. In order to estimate the location of modes and antimodes, as well as the value of the density in these points under H0, the kernel density estimator (1.1.1) with critical bandwidth ˆ fhkis a good candidate. In addition, if fis also a bounded density with bounded support [a, b]1, twice–differentiable in (a, b) and limx→a+f0(x)>0, limx→b−f0(x)<0, the result in Mammen et al. (1992) holds and the critical bandwidth is of order n−1/5, being then optimal for estimating f. Unfortunately, a suitable bandwidth for estimating fcannot be used for effectively estimating f00. In fact, the optimal bandwidth for estimating the second derivative of the density is of order n−1/9. Therefore, a simple b diobtained as the ratio |b f00 hk(bxi)|/b f3 hk(bxi), where bxiis the estimated location of the mode (antimode), will not estimate correctly di. For this reason, the kernel density estimator with critical bandwidth hkcannot be used as calibration function g, which needs to provide a good estimator of di. Under 1When [a, b] is known, the critical bandwidth proposed by Hall and York (2001) can be used. 34 CHAPTER 2. MODE TESTING FOR LINEAR DATA the conditions of Mammen et al. (1992), such an estimator can be calculated using the ratio b di=|b f00 hPI (bxi)|/b f3 hk(bxi), where hPI is a plug–in bandwidth. In this thesis, the plug–in rule for the second derivative will be obtained deriving the Asymptotic Mean Integrated Squared Error (AMISE) and replacing fin its expression using a two–step procedure (see, for example, Wand and Jones, 1995, Ch. 3). The calibration function gwill be obtained by modifying b fhkin a neighbourhood of {x:b f0 hk(x) = 0}, being such values a finite collection (see Silverman, 1981), having estimated density positive. The form of this calibration function gis given in (2.1.9) and explained below. Depending on the nature of these points, two modifications in their neighbourhood will be done. If the point bxiis a mode or an antimode of b fhkthen the modification will preserve its location and its estimated density value, and its second derivative will satisfy g00(bxi) = b f00 hPI (bxi)2. For i= 1,...,(2k−1), this will be the role of the function J, in the neighbourhood (ri,si) of bxifor guaranteeing |g00(bxi)|/g3(bxi) = b diand satisfying (RC4). The second modification will remove the tsaddle points of b fhk, denoted as ζp, with p={1, . . . , t}. In order to ensure that g satisfies (RC3), in a neighbourhood (z(2p−1), z(2p)) of ζpa function Lwill be introduced. Finally, as b fhkverifies (RC1) and (RC2), and the modifications of the functions Jand Lwill be done in bounded neighbourhoods and preserving condition (RC1), gwill satisfy also these two conditions. Then, the calibration function gwill be constructed as follows: g(x;hk, hPI,ς) =          J(x;bxi, hk, hPI, ςi) if x∈(ri,si) for some i∈ {1, . . . , (2k−1)}, L(x;z(2p−1), z(2p), hk) if x∈(z(2p−1), z(2p)) for some p∈ {1, . . . , t}, and ζp/∈(ri,si) for any i∈ {1, . . . , (2k−1)}, b fhk(x) otherwise, (2.1.9) where ςhas kcomponents ςi∈(0,1/2), with i= 1, . . . , k, determining at which height of the kernel density estimation the modification in the neighbourhood of the modes or antimodes will be done. Values of ςiclose to 0 imply a modification in an “small” neighbourhood around the mode or antimode. Note that a little abuse of notation was made as gwill depend on the function b fhk(not only on hk) and on the values b f00 hPI (bxi), for i∈ {1, . . . , (2k−1)}. An example of the effect of gcan be seen in Figure 2.3 and it will be further described later on. Before defining function J(and also function L), to ensure that ghas continuous derivative, a link function lmust 2Note that, although asymptotically the sign of b f00 hPI (bxi) is always correct, in the finite–sample case, it may not be negative in the modes or positive in the antimodes. In that case, an abuse of notation will be done, denoting as hPI to the critical or other plug–in bandwidth in order to guarantee that the sign of this second derivative remains correct. 2.1. A REVIEW ON MULTIMODALITY TESTS 35 be introduced: l(x;u, v, a0, a1, b0, b1) = a0−a1 2 1+2x−u v−u3 −3x−u v−u2!exp 2(x−u)b0 a0−a1+ +a0−a1 2 2x−u v−u3 −3x−u v−u2!exp 2(v−x)b1 a0−a1+a0+a1 2, (2.1.10) where a06=a1and v > u. Two issues must be noticed in this function. First, it will allow for a smooth connection between other two functions, being uand vthe starting and the ending points where the link function is used, a0and a1the values of the connected functions on these points and b0and b1their first derivative values. Second, if the sign of b0,b1and (a1−a2) is the same, then the first derivative of l will not be equal to 0 for any point inside the interval [u, v]. The form of the function Jis given in equation (2.1.11) and its construction guarantees that bxiis the unique point in which the derivative is equal to 0 in the neighbourhood where it is defined. The construction of Jis achieved with the K function defined below and properly linked with the link function (2.1.10) to preserve (RC1). The Kfunction is defined as follows K(x;bxi,pi,qi, ηi) = pi 1 + δix−bxi ηi2!η2 i δi·qi 2pi , being δia value indicating if bxiis a mode (δi=−1) or an antimode (δi= 1). The value ηiwill be defined later and it will depend on ςi. The second derivative of this function exists and is H¨older continuous in (bxi−ηi/2,bxi+ηi/2). The following equalities are also satisfied: K(bxi;bxi,pi,qi, ηi) = piand K00(bxi;bxi,pi,qi, ηi) = qi. Then, denoting as ρi= (bxi,b fhk(bxi),b f00 hPI (bxi)), the Jfunction can be defined as follows J(x;bxi, hk, hPI, ςi) =          lx;ri,vi,b fhk(ri),K(vi;ρi, ηi),b f0 hk(ri),K0(vi;ρi, ηi)if x∈(ri,vi), K(x;ρi, ηi)if x∈[vi,wi], lx;wi,si,K(wi;ρi, ηi),b fhk(si),K0(wi;ρi, ηi),b f0 hk(si)if x∈(wi,si), (2.1.11) being vi=bxi−ηi/2 and wi=bxi+ηi/2. As it was mentioned, the function Jdescribed in (2.1.11) (and hence also the calibration function g) depends on the constant ςi∈ (0,1/2). Ordering the modes and denoting as bx0=−∞ and bx(2k)=∞, that is −∞ =bx0<bx1< . . . < bx2k−1<bx2k=∞, the remaining unknowns values in (2.1.11) will be obtained as follows. First, it is necessary to decide at which height ϑithe modification in b fhkis done. For values of ςiclose to 0, ϑiwill be close to b fhk(bxi); while for values close to 0.5, ϑiwill be in the middle point between b fhk(bxi) 36 CHAPTER 2. MODE TESTING FOR LINEAR DATA and the highest (or lowest if bxiis an antimode) value of b fhkin the two closest modes or antimodes (bxi−1and bxi+1). Second, once the height is decided, riand siwill be the left and the right closest points to bxiat which the density estimation is equal to ϑi. Third, in order to link correctly the Kfunction, it is necessary to define ηi ensuring that K(bxi±ηi/2; ρi, ηi) will be higher (lower in the antimodes) than ϑi. With this objective, ηiis chosen in such a way that K(bxi±ηi/2; ρi, ηi) is near b fhk(bxi) and as close as possible to the middle point between ϑiand b fhk(bxi). Also, the value ηiwill ensure that the neighbourhood [vi,wi] in which Kis defined is inside (ri,si). Finally, b f0 hkmust be different to 0 in the four points (ri,vi,wiand si) where the two link functions are employed. An example of the modifications achieved by the J function in the modes and antimodes of b fhkis shown in Figure 2.3 and the complete characterization is provided below ϑi=b fhk(bxi) + δi·ςi·min |b fhk(bxi)−b fhk(bxi−1)|,|b fhk(bxi)−b fhk(bxi+1)|, ri= inf{x:x > bxi−1, δi·b fhk(x)≤δi·ϑiand b f0 hk(x)6= 0}, si= sup{x:x < bxi+1, δi·b fhk(x)≤δi·ϑiand b f0 hk(x)6= 0}, ηi= sup{γ:γ∈(0,min(bxi−ri,si−bxi)), δiK(bxi+γ/2; ρi, γ)≤δi(b fhk(bxi) + ϑi)/2 and b f0 hk(bxi±γ/2) 6= 0}. (2.1.12) In order to proceed with the modification achieved with the Lfunction, assume that this estimator has tsaddle points ζp, with p= 1, . . . , t. Define as ξ= min{|x−y|: x, y ∈(ζ1, . . . , ζt)∪(r1,s1,...,r2k−1,s2k−1)}. Then, if ζpis not inside the interval where the Jfuntions are defined, the neighbourhood used to remove the stationary and turning points will be delimited by z(2p−1) =ζp−$ξ and z(2p)=ζp+$ξ, with $∈(0,1/4). In the simulation study, the value of $will be taken close enough to 0 to avoid an impact in the value of the integral associated to g. Once these points are calculated, the saddle points can be removed from gwith the link function by taking Lequal to L(x;z(2p−1), z(2p), hk) = l(x;z(2p−1), z(2p),b fhk(z(2p−1)),b fhk(z(2p)),b f0 hk(z(2p−1)),b f0 hk(z(2p))). (2.1.13) To construct the calibration function, values of ςi∈(0,1/2), for i∈ {1, . . . , (2k− 1)}must be fixed. Then, using the Jfunction (2.1.11) with the values given in (2.1.12) and the Lfunction (2.1.13), the function gdefined in (2.1.9) satisfies the specified regularity conditions and |g00(bxi;hk, hPI,ς)|/g3(bxi;hk, hPI,ς) converges in probability to di. With this modification the calibration function also preserves the structure of 2.1. A REVIEW ON MULTIMODALITY TESTS 37 Figure 2.3: Sample of n= 1000 observations from model M16 (Appendix A). Dotted grey line: b fh2. Solid line: gfunction. Dashed line: support of J(·;bxi, h2, hPI,0.4), with i= 1,2,3, Dot–dashed line: support of K(·;ρ3, η3). Left: in the support (−0.2,1.1). Right: in a neighbourhood of the mode bx3. the data under the hypothesis that fhas kmodes, which will help to have better results in the practice. However, this calibration function gmay not be a density function since q(ς) = Z∞ −∞ g(x;hk, hPI,ς)dx (2.1.14) may not be equal to 1. To guarantee that gis indeed a density, a possible approach consists in proceeding with a search of values for ςsuch that q(ς) is equal to 1, which is guaranteed since lim ςi→0+;∀i∈(1,...,2k−1) q(ς) = Z∞ −∞ b fhk(x)dx = 1. For convenience, in the simulation study, the employed approach will be followed using ςiclose enough to 0 (∀i∈ {1,...,2k−1}) in order to avoid an impact on the integral value. To conclude this section it is important to mention the necessary condition of the bounded support in order to ensure the results of Mammen et al. (1992). Although, in general, it will not have an important effect, as shown in Section 2.2, in case of a known support, hHY,k can be used to create the calibration function. An alternative approach for known support is given below in Section 2.1.3.1. Once a conclusion 38 CHAPTER 2. MODE TESTING FOR LINEAR DATA about the number of modes is made, when both the estimation of the modes and antimodes is needed, this known support plays an special role. In Section B.1 of the Appendix, it is shown that ˆ fhHY,k provides a good estimation of their localization. It should be also reminded that the methods proposed by Silverman (1981) and Fisher and Marron (2001) can be extended from unimodality to test a general null hypothesis as H0:j≤k. Nevertheless, the proposal presented in this work just allows to test H0:j=kvs. H0:j > k. The reason why the k–modal test should not be used when the true underlying density has less than kmodes is that the test statistic in the bootstrap resamples converge in distribution to a random variable, depending only on the values b diwith i= 1,...,(2k−1) (see Cheng and Hall, 1998). When j < k, in the calibration function g, there exist (2k−2j) turning points that they will not converge to any fixed value depending on the real density function. As the (asymptotic) distribution of the test statistic in the bootstrap resamples depends also on this (2k−2j) values, one would expect that the sample distribution of the test statistic will not be correctly approximated with the bootstrap resamples. Testing H0:j=kinstead of H0:j≤kis not, in general, an important limitation for practical purposes. As it will be seen in the stamp example in Section 2.3, the usual procedure is to perform a stepwise algorithm starting with one mode and, if the null hypothesis is rejected, increasing the number of modes in the null hypothesis by one until there is no evidences for rejection. Despite this note of caution, in Section 2.2, it is showed that, generally, testing H0:j=k, when j < k, reports also good calibration results. 2.1.3.1 New proposal when the support is known When the bounded support is known, an alternative approach for the new proposal can be used in order to get better results in practice. This new proposal consists in replacing the critical bandwidth of Silverman (1981) for the one of Hall and York (2001) in the definition of the calibration function g. If the number of modes of gin the entire support is equal to kwhen testing H0:j=k(with (k−1) antimodes) in [a, b], then no more changes are needed. If modes appear outside [a, b], then the link function (2.1.10) can be used in order to preserve the required regularity conditions. Denoting as a < bx1< . . . < bx2k−1< b, being bx1and bx2k−1modes, the points bx0and bx2k, needed to obtain the values in (2.1.12), will be redefined to remove the modes outside [a, b]. If there are modes lower than bx1, then bx0= min{x:x≥aand b f0 hHY,k (x)>0}and if there are modes greater than bx2k−1, then bx2k= max{x:x≤band b f0 hHY,k (x)<0}. Once this change is done, two extra values, a<bx0and b>bx2k, are needed in order to use the link function. The steps to obtain these two values will be defined later and, 2.1. A REVIEW ON MULTIMODALITY TESTS 39 from them, the calibration function in (2.1.9) can be modified in its tails to define g(x;hHY,k, hPI,ς,a,b) as follows                              0if x≤aand b fhHY,k has modes lower than a, l(x;a,bx0,0,b fhHY,k (bx0),0,b f0 hHY,k (bx0)) if x∈(a,bx0) and b fhHY,k has modes lower than a, J(x;bxi, hHY,k, hPI, ςi) if x∈(ri,si) for some i∈ {1, . . . , (2k−1)}, L(x;ζp, hHY,k) if x∈(z(2p−1), z(2p)) for some p∈ {1, . . . , t}, and ζp/∈(ri,si) for any i∈ {1, . . . , (2k−1)}, l(x;bx2k,b,b fhHY,k (bx2k),0,b f0 hHY,k (bx2k),0) if x∈(bx2k,b) and b fhHY,k has modes greater than b, 0 if x≥band b fhHY,k has modes greater than b, b fhHY,k (x) otherwise, the functions Jand Lare defined as in Section 2.1.3, replacing the kernel density estimator b fhkby b fhHY,k and changing the values of bx0and bx2kas it was pointed out. The neighborhood in which the Jfunctions are defined is chosen by the same method as in the approach described in the main text. To guarantee that the calibration function is a density, it is also necessary to select correctly the values of a(if b fhHY,k has modes lower than a) and b(if it has modes greater than b) to obtain an integral equal to one. An option is to employ aand bsatisfying Zbx0 −∞ g(x;hHY,k, hPI,ς,a,b)dx+R∞ bx2kg(x;hHY,k, hPI,ς,a,b)dx = Zbx0 −∞ b fhk(x)dx +Z∞ bx2kb fhk(x)dx. (2.1.15) It may happen that the equality (2.1.15) is not satisfied for any pair (a,b), being a∈(−∞,bx0) and b∈(bx2k,∞). In this case, the calibration function can be divided by the normalizing constant to correct the value of the integral. Another alternative can be to take other values of bx0<bx1and bx2k>bx2k−1, such as b f0 hHY,k (x)>0, for all x∈[bx0,bx1), and b f0 hHY,k (x)<0, for all x∈(bx2k−1,bx2k]. The approach considered in the simulation study (when the support is known) is search the value of ain the interval [bx0−b+a, bx0) and the value of bin (bx2k,bx2k+b−a]. If for all the possible values of aand bthe integral q2=Z∞ −∞ g(x;hHY,k, hPI,ς,a,b)dx, is not equal to 1, then the solution is take aand bin such a way that q2is as close as possible to one and, then, employ the quotient g(·;hHY,k, hPI,ς,a,b)/q2as the calibration function. 40 CHAPTER 2. MODE TESTING FOR LINEAR DATA 2.2 Simulation study The aim of the following simulation study is to compare the level accuracy and power of the different proposals presented in Section 2.1. With this goal, samples of size n= 50, n= 200 and n= 1000 (n= 100 instead of n= 1000 in power studies) were drawn from twenty six different distributions, eleven of them unimodal (models M1–M10 and M26), ten bimodal (models M11–M20) and five trimodal (models M21–M25), see Appendix A. For each choice of sampling distribution and sample size, 500 realizations of the sample were generated. Conditionally on each of those samples, for testing purposes, 500 resamples of size nwere drawn from the population. Tables 2.1–2.7 report the percentages of rejections for significance levels α= 0.01, α= 0.05 and α= 0.10 under different scenarios: testing unimodality vs. multimodality (Table 2.1 and 2.2); testing bimodality against more than two modes when the true distribution has two modes (Table 2.4) or one (Table 2.5) mode and their power analysis (respectively Tables 2.3 and 2.6). The procedures considered include the proposals by Silverman (1981) (SI), Fisher and Marron (2001) (FM), Hall and York (2001) (HY), Hartigan and Hartigan (1985) (HH), Cheng and Hall (1998) (CH) and the new proposal (NP) in this manuscript. Note that for testing H0:j= 2, only SI, FM and NP can be compared. For the critical bandwidth test HY, the two proposed methods for computing λαhave been tried, with very similar results. The ones reported in this section correspond to a polynomial approximation for λα.I= [0,1] is used both for HY and for NP, when the interval containing the modes is assumed to be known (Tables 2.8 and 2.7). Further computational details are included in Appendix C.2. Testing unimodality vs. multimodality. From the results in Tables 2.1 and 2.2, it is concluded that SI is quite conservative, since even for high sample sizes, the percentage of rejections is below the significance level, being quite close to 0 even for α= 0.10. An exception to this conservative behaviour can be found in models M5 (with outliers in tails) and M6 (with a flat tail), finding percentages of rejections between 0.03 and 0.052 for n= 200 or n= 1000, when α= 0.10. Regarding FM, a systematic behaviour cannot be concluded: the percentage of rejections is above the significance level for models M1, M5, M7, M9 or M10, but it can be also below the true level for M2, M4 or M8. The behaviour of HY is quite good when using I= [0,1] for the different distributions and large sample sizes. For n= 1000, the percentage of rejections is quite close to α, except for model M5 (for α= 0.05, below level) and for models M6 and M7 (for α= 0.10, above level). However, the percentage of rejections is usually below the significance level for small sample sizes. Exceptions to this general pattern are found for model M1 (n= 200), M3 (n= 50) and M10 (n= 200), where percentage of rejections is close to αand models M3 (n= 200), M6 (n= 200) and M7 with percentages 2.2. SIMULATION STUDY 41 above α. Nevertheless, it should be kept in mind that the support where unimodality is tested must be known. Similarly to SI, the results obtained with HH are quite conservative. For instance, for n= 1000, even taking α= 0.10, the percentage of rejections is always below 0.002. In simple models, the CH proposal has a correct calibration, although slightly conservative in some cases, such as for model M4 (n= 1000), M5 (n= 200) and M8 (n= 200). As expected by Cheng and Hall (1998), the parametric calibration distributions do not capture, for example, the skewness and this affects the second– order properties in more complex models. This effect is reflected in the asymmetric M3, M7 and M10 (n= 1000), or model M9, where the percentage of rejections is below α, and for model M6 where is considerably higher than the significance level. Finally, regarding the new proposal NP, it seems that the calibration is quite satisfactory, even for complicated models, with a slightly conservative performance for M3 (n= 200), M4 (n= 1000) or M7 (n= 200), being this effect more clear for model M9. The unique scenario where the percentage of rejections is above αis for M6 with n= 200, but this behaviour is corrected when increasing the sample size. Although the performance is better for higher sample sizes (n= 200 or n= 1000), in some cases, such as M9 or M16 (in the bimodal case), even for n= 1000, a percentage of rejections close to αis hard to get. In this difficult cases, the knowledge of the support can be used for obtaining better results as it was reported in Table 2.8, where the percentage of rejections is close to αfor the sample sizes n= 200 and n= 1000. Regarding power behaviour (and just commenting on the three methods which exhibit a correct calibration: HY, CH and NP), results are reported in Table 2.3: none of the proposals is clearly more powerful. For instance, for model M13, HY clearly detects (even for n= 50) the appearance of the second small mode, whereas the other approaches do not succeed in doing so. For models M11, M12, M14 and M15 (n= 50), CH presents the highest empirical power and HY shows the lowest one. Assessing bimodality. For testing H0:j= 2, Hall and York (2001) prove that, even knowing the density support, SI cannot be consistently calibrated by a bootstrap procedure, similar to the one used for the unimodality vs. multimodality test. Hence, in some models such as M15 or M16 (n= 50), the percentage of rejections is close to α, whereas the conservative behaviour observed in the unimodality test, is just perceived in M14 and M19, for large sample sizes. Also, in this case, there is a model where the percentage of rejections is considerably higher than the significance level, M20 (n= 200 and n= 1000), being this bimodal model similar to the conservative M19, just generating some atypical data in the tails. FM presents again an erratic behaviour: for M17 (except for n= 200), M18 or M19, the percentage of rejections 48 CHAPTER 2. MODE TESTING FOR LINEAR DATA α0.01 0.05 0.10 0.01 0.05 0.10 0.01 0.05 0.10 M21 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) SI FM NP n=50 0(0) 0.026(0.014) 0.058(0.020) 0.010(0.009) 0.052(0.019) 0.112(0.028) 0.014(0.010) 0.054(0.020) 0.112(0.028) n=100 0.002(0.004) 0.056(0.020) 0.154(0.032) 0.044(0.018) 0.148(0.031) 0.228(0.037) 0.020(0.012) 0.092(0.025) 0.168(0.033) n=200 0.024(0.013) 0.104(0.027) 0.208(0.036) 0.076(0.023) 0.202(0.035) 0.278(0.039) 0.050(0.019) 0.134(0.030) 0.194(0.035) M22 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) SI FM NP n=50 0.002(0.004) 0.006(0.007) 0.022(0.013) 0.036(0.016) 0.142(0.031) 0.300(0.040) 0.080(0.024) 0.256(0.038) 0.402(0.043) n= 100 0(0) 0.026(0.014) 0.180(0.034) 0.248(0.038) 0.610(0.043) 0.804(0.035) 0.266(0.039) 0.542(0.044) 0.706(0.040) n=200 0.052(0.019) 0.500(0.044) 0.806(0.035) 0.740(0.038) 0.960(0.017) 0.982(0.012) 0.890(0.027) 0.956(0.018) 0.980(0.012) M23 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) SI FM NP n=50 0(0) 0.068(0.022) 0.206(0.035) 0.054(0.020) 0.196(0.035) 0.338(0.041) 0.050(0.019) 0.204(0.035) 0.336(0.041) n=100 0.082(0.024) 0.486(0.044) 0.684(0.041) 0.334(0.041) 0.658(0.042) 0.780(0.036) 0.326(0.041) 0.624(0.042) 0.746(0.038) n=200 0.570(0.043) 0.780(0.036) 0.840(0.032) 0.832(0.033) 0.906(0.026) 0.938(0.021) 0.722(0.039) 0.878(0.029) 0.934(0.022) M24 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) SI FM NP n=50 0(0) 0.008(0.008) 0.038(0.017) 0.004(0.006) 0.030(0.015) 0.082(0.024) 0.008(0.008) 0.060(0.021) 0.134(0.030) n=100 0(0) 0.006(0.007) 0.026(0.014) 0.012(0.010) 0.052(0.019) 0.108(0.027) 0.034(0.016) 0.110(0.027) 0.176(0.033) n=200 0.006(0.007) 0.040(0.017) 0.108(0.027) 0.050(0.019) 0.160(0.032) 0.288(0.040) 0.096(0.026) 0.232(0.037) 0.334(0.041) M25 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) SI FM NP n=50 0.002(0.004) 0.046(0.018) 0.138(0.030) 0.068(0.022) 0.200(0.035) 0.312(0.041) 0.018(0.012) 0.098(0.026) 0.148(0.031) n=100 0.004(0.006) 0.064(0.021) 0.176(0.033) 0.170(0.033) 0.344(0.042) 0.462(0.044) 0.098(0.026) 0.232(0.037) 0.320(0.041) n=200 0.008(0.008) 0.084(0.024) 0.198(0.035) 0.190(0.034) 0.404(0.043) 0.552(0.044) 0.106(0.027) 0.248(0.038) 0.356(0.042) Table 2.6: Percentages of rejections for testing bimodality, with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. α0.01 0.05 0.10 0.01 0.05 0.10 0.01 0.05 0.10 M26 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) NP (unknown support) NP (unknown support) NP (known support) n=50 0.008(0.008) 0.032(0.015) 0.066(0.022) 0.004(0.006) 0.028(0.014) 0.048(0.019) 0.002(0.004) 0.022(0.013) 0.076(0.023) n=200 0(0) 0.022(0.013) 0.068(0.022) 0.006(0.007) 0.030(0.015) 0.074(0.023) 0.004(0.006) 0.034(0.016) 0.062(0.021) n=1000 0.004(0.006) 0.040(0.017) 0.080(0.024) 0(0) 0.020(0.012) 0.054(0.020) 0(0) 0.026(0.014) 0.052(0.019) Table 2.7: Percentages of rejections for testing unimodality (first column) and bimodality (second and third column), with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 2.3. THE 1872 HIDALGO STAMP ISSUE OF MEXICO 49 α0.01 0.05 0.10 NP (known support) n= 50 M9 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 x f(x) 0(0) 0.014(0.010) 0.034(0.016) n= 200 0.012(0.010) 0.052(0.019) 0.090(0.025) n= 1000 0.010(0.009) 0.040(0.017) 0.092(0.025) n= 50 M16 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 2.5 x f(x) 0.002(0.004) 0.024(0.013) 0.076(0.023) n= 200 0.008(0.008) 0.058(0.020) 0.126(0.029) n= 1000 0.016(0.011) 0.038(0.017) 0.080(0.024) Table 2.8: Percentages of rejections for testing unimodality (in model M9) and bimodality (in model M16), with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 2.3 The 1872 Hidalgo stamp issue of Mexico As explained in Section 1.3.1, the value of stamps depends on its scarcity, and thickness is determinant in this sense. However, in general, the designation of thick, medium or thin stamps is relative and can only refer to a particular stamp issue. Otherwise, making uniform categories for all stamp issues may lead to inaccurate classifications. In addition, there is not such a differentiation between groups available in stamps catalogues, leaving this classification to a personal subjective judgement. The importance of establishing an objective criterion specially appears in stamp issues printed on a mixture of paper types, with possible differences in their thickness. A stamp issue where the problem of determining the number of different groups of stamps appears is in the 1872 Hidalgo issue. First, for this particular issue, and in general in the Mexican ones, it is known that the handmade paper presents a lot of variability in the thickness of the paper. Second, since of scarcity of ordinary white wove paper, other types of paper were used to produce the Hidalgo issue. A small quantity of “vertically laid” paper, a fiscal type of white wove paper denominated Papel Sellado (some of them were watermarked vertically), other type of white wove, the La Croix–Freres of France (some of them with a watermark of LA+-F) and also another unwatermarked white wove paper might also have been used. It is estimated than the watermark of Papel Sellado can appear in between 6 and 18 stamps in each sheet of 100. For the La Croix–Freres watermark, it is estimated that the symbol appears in between 4 and 10 stamps if the sheet was watermarked, and some authors suggested that this watermark appears only once of every 4 sheets. In order to get more information about this particular problem and to obtain some further references, see Izenman and Sommer (1988). 50 CHAPTER 2. MODE TESTING FOR LINEAR DATA This particular example has been explored in several references in the literature for determining the number of different groups of stamps. From a nonparametric point of view, some examples of its utilization for mode testing can be shown in Efron and Tibshirani (1994, Ch. 16), Izenman and Sommer (1988) or in Fisher and Marron (2001). Also it was analysed using nonparametric exploratory tools in Wilson (1983), Minnotte and Scott (1993) and in Chaudhuri and Marron (1999). Some parametric studies of the 1872 Hidalgo issue can be found in Basford et al. (1997) or in McLachlan and Peel (2000, Ch. 6). Taking a subsample of 437 stamps on white wove, Wilson (1983) represented a histogram and the conclusion was that only two kinds of paper were used, the Papel Sellado and the La Croix–Freres, and that there was not a third kind of paper. Izenman and Sommer (1988) revisited the example considering a more complete collection, with 485 stamps. A histogram with the same parameters as those used by Wilson (1983) (same starting point and bin width) is shown in Figure 2.1 (top–left panel), revealing the same features as those noticed in the original reference. Two groups are also shown by a kernel density estimator, shown in the same plot, considering a Gaussian kernel and a rule of thumb bandwidth (see Wand and Jones, 1995, Section 3.2). However, both approximations (histogram and kernel density estimator) depend heavily on the bin width and bandwidth, respectively. Specifically, the use of an automatic rule for selecting the bandwidth value (focused on the global estimation of the entire density function) does not guarantee an appropriate recovery of the modes. In fact, using another automatic rule as the plug–in bandwidth (Figure 2.1, bottom–left panel), nine modes are observed. A histogram with a smaller bin width is also included in this plot, exhibiting apparently more modes than the initial one. Given that the exploratory tools did not provide a formal way of determining if there are more than two groups, Papel Sellado and La Croix–Freres, Izenman and Sommer (1988) decided to employ the multimodality test of Silverman (1981). Note that for this purpose, just FM, SI and the new proposal NP can be employed, and the first two proposals present a poor calibration, as shown in the simulation study. Results from Izenman and Sommer (1988), applying SI with B= 100, are shown in Table 2.9 (note that with B= 500, different p–values are obtained). For α= 0.05, the conclusions are the same, except in the crucial case of testing H0:j≤2, where for B= 500, there are no evidences to reject the null hypothesis. These differences may be caused by the approximations implemented by Izenman and Sommer (1988) to obtain the critical bandwidth. Both Efron and Tibshirani (1994, Ch. 16) (using B= 500 bootstrap replicates) and Salgado-Ugarte et al. (1998) (employing B= 600) obtained similar results to ours. Hence, the null hypothesis must not be rejected when the hypothesis is that the distribution has at most two modes, but it has to be rejected when H0is that the distribution has at most six modes. This strange behaviour also 2.3. THE 1872 HIDALGO STAMP ISSUE OF MEXICO 51 k123456789 SI B=100 00.04 0.06 0.01 0 0 0.44 0.31 0.82 B=500 0.006 0.348 0.092 0.016 0.006 0.002 0.494 0.308 0.616 FM 0 0.04 0 0 0 0 0.06 0.01 0.06 NP 0 0.022 0.004 0.506 0.574 0.566 0.376 0.886 0.808 Table 2.9: P–values obtained using different proposals for testing k–modality, with kbetween 1 and 9. Methods: SI, FM and NP. For SI method, B= 100 (first row; Izenman and Sommer, 1988) and B= 500. happens in Izenman and Sommer (1988) analysis, when testing H0:j≤3 and H0:j≤6. Izenman and Sommer (1988) suggested non–rejecting the null hypothesis the first time that the p–value is higher than 0.4. The consideration of a flexible rule for rejecting the null hypothesis is justified by the fluctuations in the p–values of SI and, as Izenman and Sommer (1988) mentioned, by the “conservative” nature of this test. Under this premise, the result when applying SI would be that the null hypothesis is rejected until it is tested H0:j≤7. Hence, Izenman and Sommer (1988) conclude that the number of groups in the 1872 Hidalgo Issue is seven. As shown in Section 2.2, SI does not present a good calibration and sometimes it can be also anticonservative. It is not surprising that SI behaves differently when testing H0:j≤kfor k= 2,3, with respect to the rest of cases until k= 7. Since NP has a good calibration behaviour, even with “small” sample sizes, this method is going to be used, first for testing the important case H0:j= 2 vs. Ha:j > 2 and then to figure out how many groups are there in the 1872 Hidalgo Issue. The computation of the excess mass statistic requires a non–discrete sample and the original data (denoted as X) contained repeated values, the artificial sample Y=X+Ewill be employed for testing the number of modes, where Eis a sample of size 485 from the U(−5·10−4,5·10−4) distribution. This modification of the data was also considered by Fisher and Marron (2001). The p–values obtained in their studio (using B= 200 bootstrap replicates) are shown in Table 2.9: it is not clear which conclusion has to be made. They mentioned that their results are consistent with the previous studies, detecting 7 modes. But it should be noticed that, as shown in the simulation, FM does not present a good calibration behaviour. Finally, the p–values obtained with NP are also shown in Table 2.9, with B= 500. Similar results can be obtained employing the interval I= [0.04,0.15] in NP with 52 CHAPTER 2. MODE TESTING FOR LINEAR DATA known support, as Izenman and Sommer (1988) notice that the thickness of the stamps is always in this interval I. Employing a significance level α= 0.05 for testing H0:j= 2, the null hypothesis is rejected. This null hypothesis is rejected until k= 4, and then there is no evidences to reject H0employing greater values of k. Then, applying our new procedure, the conclusion is that the number of groups in the 1872 Hidalgo Issue is four. In order to compare the results obtained by Izenman and Sommer (1988) and the ones derived applying the new proposal, two kernel density estimators, with Gaussian kernel and critical bandwidths h4and h7are depicted in Figure 2.1 (bottom–right panel). Izenman and Sommer (1988) conclude that seven modes were present, and argued that the stamps could be divided in, first, three groups (pelure paper with mode at 0.072 mm, related with the forged stamps; the medium paper in the point 0.080 mm; and the thick paper at 0.090 mm). Given the efforts made in the new issue in 1872 to avoid forged stamps, it seems quite reasonable to assume that the group associated with the pelure paper had disappeared in this new issue. In that case, the asymmetry in the first mode using h4can be attributed to the modifications in the paper made by the manufacturers. Also, this first and asymmetric group, justifies the application of nonparametric techniques to determine the number of groups. In the Section 7 of Izenman and Sommer (1988) and in other references using mixtures of Gaussian densities to model this data (see, for example, McLachlan and Peel, 2000, Ch. 6), it is shown that these parametric techniques have problems in capturing this asymmetry, and they always determine that there are two modes in this first part of the density, one near the point 0.07 mm and another one near 0.08 mm. For the two modes near the points 0.10 and 0.11 mm, both corresponding to stamps produced in 1872. As Izenman and Sommer (1988) noticed, it seems that the stamps of 1872 were printed on two different paper types, one with the same characteristics as the unwatermarked white wove paper used in the 1868 issue, and a second much thicker paper that disappeared completely by the end of 1872. Using this explanation, it seems quite reasonable to think that the two final modes using h4, corresponds with the medium paper and the thick paper in this second block of stamps produced in 1872. Finally, for the two minor modes appearing near 0.12 and 0.13 mm, when h7is used, Izenman and Sommer (1988) do not find an explanation and they mention that probably they could be artefacts of the estimation procedure. This seems to confirm the conclusions obtained with our new procedure. The reason of determining more groups than the four obtained with our proposal, seems to be quite similar to that of the model M20 in the simulation study (Section 2.2). This possible explanation is that the spurious data in the right tail of the last mode are causing the rejection of H0, when SI is used. 2.4. SEROPOSITIVITY IN MALARIA ERADICATION 53 2.4 Seropositivity in malaria eradication As mentioned in the Introduction, multimodality tests are not only important per se, but also as a previous step for applying other procedures. An example where this happens is found when one wants to model and determine a threshold dividing the healthy people and those ones infected with a malaria disease (see Section 1.3.2). Classical methods to understand the degree of malaria exposure in a population (as, for example, catch the mosquitoes when they attempt to land on the exposed limbs of field workers), can be quite laborious, specially in low transmission populations. For that reason, in the last few years, alternative indicators such as those ones based on antibodies against the malaria antigens are being used. The problem of using such data is that there is not a fixed rule determining from which quantity on an individual is infected (seropositive) or not (seronegative). Even more, in a given infected population, there can be different kinds of seropositive subpopulations. For that reason, in each population it is necessary to analyse each antigen in order to understand the true malaria exposure. Identifying seropositivity is simplified in studies where there is a training set of individuals with known serological status. In such settings, one can use this information to define a cutoff for seropositivity, e.g., using the sample mean plus three times the estimated standard deviation of the seronegative population (Drakeley et al., 2005) or constructing a mixture model from this prior information (Vounatsou et al., 1998). The problem is that, in general, such training set is not easily available in endemic areas and, therefore, the serological classification of the individuals collected from those areas is more ambiguous. When nothing about the serological status of the individuals is known, a first step, to understand the degree of malaria exposure in a given population, is analyse if there are at least two subpopulations (one seropositive and at least another seronegative). For seeing that, one can test if the number of modes in an antibody sample is equal or greater than one. If the number of modes is equal to one, two things may be happen, that the antibodies against the malaria antigens are similar in both subpopulations or that the malaria is (almost) eliminated from the population. If more than one mode appears, it means that at least one seropositive group is still in the population, as it is expected that the antigens of the seronegative individuals behave in a similar way and show just one peak. In this example, the multimodality tests can be employed to analyse if the malaria disease is still in the population. Then, these tools can be used as a preliminary step to adjust more complicated parametric models in order to catch the (at least) two subpopulation. Also, when the support of the antibody data is known, and a conclusion about the number of modes is obtained, due the results reported in Section B.1, one can use the kernel density estimator (1.1.2), with the critical bandwidth, to estimate the 54 CHAPTER 2. MODE TESTING FOR LINEAR DATA cutoff between the seronegative and the seropositive subpopulations. A simple approach for doing so is taking the antimode in this density estimation. In that case, the first and principal group is expected to be related with the seronegative individuals and the following groups with the seropositive ones. With illustrative purposes, the new method for testing multimodality will be applied to two case studies from malaria elimination settings: Aneityum island in Vanuatu archipelago and Chabahar in southeastern Iran. In these two regions, it is known that the Plasmodium falciparum (Pf) and Plasmodium vivax (Pv) parasites, which cause malaria diesase, co–exist. The Aneityum dataset (n= 517) is part of a larger cross–sectional survey conducted in 2009 in Tafea province, Vanuatu (Cook et al., 2010). This small island has approximately 800 permanent residents and has been under strong malaria control since 1991 (Chan et al., 2017). According to annual surveillance surveys routinely conducted in the island, there are no reports of malaria infections between 1992 and 2009 with the exception of an Plasmodium vivax occurred in the year of 2002. The Iran dataset comprises a total number of 1479 people of all ages collected from the Chabahar city and surrounding villages in the southeastern part of the country (Zakeri et al., 2016). In both datasets, antibody titers against the antibody Merozoite Surface Protein–1 for Pf or Pv were analysed. After applying the new procedure for testing unimodality (with B= 1000 bootstrap resamples), for the two antigens (PfMSP1 and PvMSP1) and in the two populations (Aneityum and Chabahar), there are no evidences (using a significance level of α= 0.10) for rejecting the null hypothesis. The obtained p–values are provided in Table 2.10. Since the null hypothesis of unimodality is not rejected in any case, in Figure 2.4, the kernel density estimation with the Silverman (1981) critical bandwidth for one mode is shown. Observing the smoothed density estimation in the four cases, one may ask if the outlier data far from the zero are altering the results of the test. For that reason, the data study was repeated assuming that the compact support was known. As expected with the results reported in Section 2.2, analysing Table 2.10 is clear that the non–rejection of the null hypothesis is not caused by this factor. Taking this compact support, the kernel density estimation with the critical bandwidth of Hall and York (2001) seems more accurate. The non–rejection of unimodality seems to be caused by the almost complete elimination of the malaria disease in both datasets and in case that it remains seropositive individuals, they are atypical data in the population. Also, observing the shape of the data, one can think about applying a transformation of the data, such us the Box–Cox method. Although in that case, one should be careful with the extracted conclusions, as there is no guarantee that more than one mode in the transformed data imply one than more group in the original data. In any case, applying this transformation, the results are consistent, as the null hypothesis is still not rejected for a significance level of α= 0.10 in all the datasets. 2.4. SEROPOSITIVITY IN MALARIA ERADICATION 55 Finally, for future analysis in other datasets, if one knows that there are still two groups in the population, a natural question is: where is the threshold dividing the two subpopulations (seronegative and seropositive)? As aforementioned, in that case, the critical bandwidth can be used for that purpose, but one must be careful with the derived conclusions as the bounded support is unknown. Thus, observing Table 2.11, in the case that the critical bandwidth of Silverman (1981) is used, in the four cases, just the highest observation is detected as seropositive. If previous information about the typical values of an antigen in the population are provided, then they can be used for a more accurate study determining the cutoff point and therefore which individuals can be seropositive, as shown in Table 2.11. 56 CHAPTER 2. MODE TESTING FOR LINEAR DATA Aneityum PfMSP1 Support R[0,100] [0,200] P–value 0.914 0.874 0.881 PvMSP1 Support R[0,100] [0,200] P–value 0.902 0.831 0.871 Chabahar PfMSP1 Support R[0,200] [0,500] P–value 0.884 0.866 0.865 PvMSP1 Support R[0,200] [0,300] P–value 0.504 0.394 0.407 Table 2.10: P–values obtained using the new proposal for testing unimodality in the antigens data (B= 1000). In the two last columns, since the support is given, the proposal provided in Section 2.1.3.1 is employed. Unimodality Bimodality Mode Mode 1 Antimode Mode 2 Aneityum PfMSP1 Support R[0,100] [0,200] R[0,100] [0,200] R[0,100] [0,200] R[0,100] [0,200] Location 28.93 21.42 24.46 24.46 20.90 23.86 225.48 90.10 134.95 283.08 99.32 147.29 PvMSP1 Support R[0,100] [0,200] R[0,100] [0,200] R[0,100] [0,200] R[0,100] [0,200] Location 32.01 23.46 27.33 31.64 22.45 24.85 401.40 83.94 136.49 443.75 89.18 158.86 Chabahar PfMSP1 Support R[0,200] [0,500] R[0,200] [0,500] R[0,200] [0,500] R[0,200] [0,500] Location 47.20 31.84 29.76 45.91 29.91 34.46 1505.30 168.18 218.52 2488.26 173.18 277.27 PvMSP1 Support R[0,200] [0,300] R[0,200] [0,300] R[0,200] [0,300] R[0,200] [0,300] Location 38.27 31.13 31.15 34.83 30.93 31.15 673.50 182.20 221.38 797.74 192.72 262.77 Table 2.11: Estimated location of the modes and antimodes in the antigens data. 2.4. SEROPOSITIVITY IN MALARIA ERADICATION 57 0 50 100 150 200 250 300 0.000 0.010 0.020 0.030 Antibody titres 0 50 100 150 0.000 0.010 0.020 0.030 Antibody titres 0 100 200 300 400 0.000 0.005 0.010 0.015 0.020 0.025 0.030 Antibody titres 0 50 100 150 0.000 0.005 0.010 0.015 0.020 0.025 0.030 Antibody titres 0 500 1000 1500 2000 2500 0.000 0.005 0.010 0.015 0.020 Antibody titres 0 50 100 150 200 250 300 0.000 0.005 0.010 0.015 0.020 Antibody titres 0 200 400 600 800 0.000 0.005 0.010 0.015 0.020 Antibody titres 0 50 100 150 200 0.000 0.005 0.010 0.015 0.020 Antibody titres Figure 2.4: Kernel density estimation with Gaussian kernel and critical bandwidth for one mode. Top: Aneityum dataset. Bottom: Chabalar observations. Left: antigen PfMSP1 (first column: in the entire range, second column: in a subinterval). Right: antigen PvMSP1 (third column: in the entire range, forth column: in a subinterval). Dotted lines: Silverman (1981) critical bandwidth. Dashed and solid lines: Hall and York (2001) critical bandwidth. Dashed lines: in the support [0,200], Aneityum datasets; [0,500], PfMSP1 in the Chabalar case; [0,300], PvMSP1 in Chabalar. solid lines: in the support [0,100], Aneityum datasets; [0,200], Chabalar samples. The different samples are represented at the bottom of the graphics with ticks. 64 CHAPTER 3. MODE TESTING IN CIRCULAR DATA Figure 3.1: Sample of n= 200 observations obtained from model MC7 (described in Appendix A). Dotted grey line: ˆ fν1. Solid line: gfunction. Dashed line: neighbourhood where the functions J(·;b θi, ν1, νPI,0.25) are defined, with i= 1,2. Dot–dashed line: neighbourhood where K(·;ρ1, η1) is defined. Left: in the support [0,2π). Right: in a neighbourhood of the mode b θ1. the generated samples, for a significance level α, the null hypothesis will be rejected if P(∆∗ n,k+1 ≤∆n,k+1|X)≥1−α. Since, in the circular case, an analogous of the critical bandwidth in the bounded support [0,2π) is obtained, the critical concentration defined in (3.2.1), one can think in adapting the ideas of Hall and York (2001) for testing unimodality. In this case, even for the most simple case, H0:j= 1, observing the simulation results it seems that the calibration of the test statistic is not independent of the underlying distribution of the data, which make sense as in this case the number of turning points of fin [0,2π) is always greater than one and, in this situation, even in the linear case (see Section 2.1.1), the test statistic cannot be directly calibrated. 3.3 Simulation study The aim of the following simulation study is to analize the performance of the bootstrap procedure to calibrate the circular excess mass statistic using the calibration 3.3. SIMULATION STUDY 65 function gdescribed in (3.2.3) for obtaining the bootstrap resamples. In this simulation study, it is also compared the level accuracy and power of our new method with the other proposal for testing multimodality for circular data, the one introduced by Fisher and Marron (2001). They propose using the U2statistic of Watson (1961) as a test statistic, that is U2=nZ2π 0Fn(x)−F0(x)−Z2π 0 (Fn(y)−F0(y))dF0(y)2 dF0(x), estimating F0(circular distribution function) employing a kernel distribution estimation with kmodes, when the problem of interest is testing H0:j≤k. In this simulation study, to estimate F0, it is used the distribution function associated to ˆ fνk and its associated distribution is used to generate the bootstrap resamples to calibrate the test statistic. With these objectives, samples of size n= 50, n= 200 and n= 1000 (n= 100 instead of n= 1000 in power studies) were drawn from 25 different distribution, ten of them unimodal (MC1–MC10), ten bimodal (MC11–MC20) and five trimodal (MC21–MC25), including unimodal symmetric models, mixtures of them and reflective asymmetric models. Those distributions are described in Appendix A. For each model (MC1–MC25) and sample size, 500 sample realizations were generated. Conditionally on each sample, for each test, 500 resamples of size nwere drawn using the calibration function of that specific test. Tables 3.1–3.6 report estimates of the nominal levels α= 0.01, α= 0.05 and α= 0.10. Tables 3.1 and 3.2 will be used to study the calibration of these tests when the problem of interest is testing H0:j= 1. Table 3.3 will be employed to analyse the power of the proposals for testing unimodality. Tables 3.4 and 3.5 will be used to study the calibration of these tests when the problem of interest is H0:j= 2. Table 3.6 will be employed to analyse the power of the proposals for testing bimodality. Analysing Tables 3.1, 3.2, 3.4 and 3.5, the poor calibration of the Fisher and Marron (2001) proposal can be observed. Even for sample size equal to 1000, sometimes, the percentage of rejections is under the significance level, as in the distributions where unimodality is tested: MC1, MC2, MC4, MC8, MC9 or MC10; or the models in which bimodality is tested: MC12, MC14, MC16 or MC17. It should be also noticed that for models MC3 and MC5 (unimodality) and MC11 and MC20 (bimodality), the percentage of rejections is above α. Studying the calibration of our new proposal (Tables 3.1, 3.2, 3.4 and 3.5), it provides, in general, a good level accuracy, except for model MC3. Even for small sample sizes (n= 50), when the null hypothesis of unimodality is tested, the percentage of rejections is close to the significance level α, except in the commented case 66 CHAPTER 3. MODE TESTING IN CIRCULAR DATA MC3, and also on models: MC1 (n= 50), MC4 (n= 200), MC8 (n= 50), MC9 (n= 200) and MC10 (n= 50), where the percentage of rejections is slightly below the significance level. In the case of testing bimodality, taking a sample size equal at least to n= 200, our proposal seems to calibrate correctly, except for model MC11 where, the percentage of rejections is slightly below α. In the case of the model MC3, the conservative behaviour observed, even for sample size equal to 1000, is corrected when considering a larger sample size (n= 2000). In Tables 3.3 and 3.6, power results are presented. Our new proposal, which seems to be the only one which is well calibrated, appears to have also good power, in terms that the percentage of rejections increases with the sample size. This method detects the clearly rejection of the null hypothesis on the bimodal model MC11 and on the trimodal models MC21 and MC22. This new proposal, also detects the small blips, for example on models MC14, MC15 (bimodal) and MC25 (trimodal), although, in other distributions, it has some troubles, when this small blips represents a low percentage of the data and also the sample size is small, like, for example, in model MC10. In the difficult cases, with almost overlapping peaks, such as on models MC13 (bimodal), MC23 (trimodal) and MC24 (trimodal), when the sample size is small, our method has some difficulties to detect the rejection of unimodality (with n= 50 in MC13) and the rejection of bimodality (with n= 100 in MC23 and with n= 50 in MC24), but the percentage of rejections increases with n. 3.4 Data analysis: detection of fire seasonality As explained in Section 3.1, determining the number of fire seasons in the different regions of the world has a special importance when the relationships between the climate and the land managment are studied. Using the data described in Section 3.4.1, the goal is to assess, in the study area, if there is one ore more fire seasons. If just one cell is considered, this problem can be tackled employing the new procedure introduced in Section 3.2, as the simulation study in Section 3.3 supports, in the finite–sample case, that the proposal presents a correct behaviour in terms of calibration and power. Since the goal is to analyse the number of fire seasons all across the world, the globe can be divided in a grid and then the proposed procedure can be applied systematically in each cell. As mentioned before, a FDR procedure is required in order to control the incorrect rejections of the null hypothesis. For performing such a correction, the spatial correlation between the p–values computed at different cells must be considered. These two issues will be solved in Section 3.4.2 and the obtained results are provided in Section 3.4.3. 3.4. DATA ANALYSIS: DETECTION OF FIRE SEASONALITY 67 α0.01 0.05 0.10 MC1 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.008(0.008) 0.028(0.014) 0.064(0.021) n=200 0.006(0.007) 0.030(0.015) 0.064(0.021) n=1000 0.004(0.006) 0.032(0.015) 0.058(0.020) Excess mass n=50 0.002(0.004) 0.022(0.013) 0.070(0.022) n=200 0.004(0.006) 0.034(0.016) 0.074(0.023) n=1000 0.006(0.007) 0.052(0.019) 0.088(0.025) MC2 0.0 0.2 0.4 0.6 0.8 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0(0) 0.022(0.013) 0.068(0.022) n=200 0(0) 0.012(0.010) 0.058(0.020) n=1000 0.004(0.006) 0.008(0.008) 0.032(0.015) Excess mass n=50 0.012(0.010) 0.054(0.020) 0.094(0.026) n=200 0.008(0.008) 0.038(0.017) 0.092(0.025) n=1000 0.006(0.007) 0.040(0.017) 0.086(0.025) MC3 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.654(0.042) 0.828(0.033) 0.886(0.028) n=200 0.708(0.040) 0.822(0.034) 0.884(0.028) n=1000 0.574(0.043) 0.718(0.039) 0.770(0.037) Excess mass n=50 0(0) 0.018(0.012) 0.038(0.017) n=200 0.002(0.004) 0.012(0.010) 0.024(0.013) n=1000 0(0) 0.014(0.010) 0.022(0.013) MC4 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0(0) 0.024(0.013) 0.066(0.022) n=200 0(0) 0.020(0.012) 0.042(0.018) n=1000 0.004(0.006) 0.024(0.013) 0.058(0.020) Excess mass n=50 0.006(0.007) 0.044(0.018) 0.114(0.028) n=200 0.008(0.008) 0.026(0.014) 0.068(0.022) n=1000 0.01(0.008) 0.046(0.018) 0.082(0.024) MC5 0.0 0.2 0.4 0.6 0.8 1.0 1.2 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.370(0.042) 0.520(0.044) 0.598(0.043) n=200 0.648(0.042) 0.776(0.037) 0.848(0.031) n=1000 0.498(0.044) 0.648(0.042) 0.728(0.039) Excess mass n=50 0.008(0.008) 0.044(0.018) 0.104(0.027) n=200 0.014(0.010) 0.036(0.016) 0.076(0.023) n=1000 0.012(0.010) 0.040(0.017) 0.074(0.023) Table 3.1: Percentages of rejections for testing unimodality, with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 68 CHAPTER 3. MODE TESTING IN CIRCULAR DATA α0.01 0.05 0.10 MC6 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.006(0.007) 0.040(0.017) 0.064(0.021) n=200 0.004(0.006) 0.044(0.018) 0.108(0.027) n=1000 0.016(0.011) 0.070(0.022) 0.132(0.030) Excess mass n=50 0.002(0.004) 0.046(0.018) 0.094(0.026) n=200 0.020(0.012) 0.058(0.020) 0.118(0.028) n=1000 0.006(0.007) 0.048(0.019) 0.098(0.026) MC7 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.002(0.004) 0.020(0.012) 0.052(0.019) n=200 0.012(0.010) 0.048(0.019) 0.096(0.026) n=1000 0.010(0.009) 0.050(0.019) 0.092(0.025) Excess mass n=50 0.008(0.008) 0.034(0.016) 0.086(0.025) n=200 0.012(0.010) 0.064(0.021) 0.114(0.028) n=1000 0.004(0.006) 0.040(0.017) 0.094(0.026) MC8 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.004(0.006) 0.030(0.015) 0.062(0.021) n=200 0.004(0.006) 0.032(0.015) 0.052(0.019) n=1000 0.002(0.004) 0.018(0.012) 0.048(0.019) Excess mass n=50 0(0) 0.024(0.013) 0.044(0.018) n=200 0.006(0.007) 0.044(0.018) 0.100(0.026) n=1000 0.010(0.009) 0.048(0.019) 0.100(0.026) MC9 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0 π 2 π 3π 2 + U2statistic n=50 0.010(0.009) 0.054(0.020) 0.114(0.028) n=200 0.008(0.008) 0.020(0.012) 0.076(0.023) n=1000 0(0) 0.010(0.009) 0.044(0.018) Excess mass n=50 0.004(0.006) 0.038(0.017) 0.064(0.021) n=200 0.008(0.008) 0.030(0.015) 0.062(0.021) n=1000 0.004(0.006) 0.036(0.016) 0.084(0.024) MC10 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0 π 2 π 3π 2 + U2statistic n=50 0.014(0.010) 0.056(0.020) 0.106(0.027) n=200 0.012(0.010) 0.064(0.021) 0.122(0.029) n=1000 0.004(0.006) 0.026(0.014) 0.064(0.021) Excess mass n=50 0(0) 0.030(0.015) 0.070(0.022) n=200 0.012(0.010) 0.038(0.017) 0.082(0.024) n=1000 0.002(0.004) 0.038(0.017) 0.086(0.024) Table 3.2: Percentages of rejections for testing unimodality, with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 3.4. DATA ANALYSIS: DETECTION OF FIRE SEASONALITY 69 α0.01 0.05 0.10 MC11 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.782(0.036) 0.920(0.024) 0.958(0.018) n=100 0.978(0.013) 0.996(0.006) 0.998(0.004) n=200 1(0) 1(0) 1(0) Excess mass n=50 0.534(0.044) 0.758(0.038) 0.854(0.031) n=100 0.914(0.025) 0.968(0.015) 0.984(0.011) n=200 0.996(0.006) 1(0) 1(0) MC12 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.338(0.041) 0.548(0.044) 0.654(0.042) n=100 0.594(0.043) 0.758(0.038) 0.830(0.033) n=200 0.880(0.028) 0.940(0.021) 0.956(0.018) Excess mass n=50 0.002(0.004) 0.040(0.017) 0.076(0.023) n=100 0.010(0.009) 0.058(0.020) 0.112(0.028) n=200 0.040(0.017) 0.116(0.028) 0.238(0.037) MC13 0.00 0.05 0.10 0.15 0.20 0.25 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.014(0.010) 0.070(0.022) 0.130(0.029) n=100 0.038(0.017) 0.124(0.029) 0.196(0.035) n=200 0.066(0.022) 0.170(0.033) 0.246(0.038) Excess mass n=50 0.022(0.013) 0.072(0.023) 0.140(0.030) n=100 0.032(0.015) 0.084(0.024) 0.156(0.032) n=200 0.044(0.018) 0.102(0.027) 0.204(0.035) MC14 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.476(0.044) 0.730(0.039) 0.840(0.032) n=100 0.828(0.033) 0.952(0.019) 0.976(0.013) n=200 0.990(0.009) 0.998(0.004) 1(0) Excess mass n=50 0.044(0.018) 0.164(0.032) 0.284(0.040) n=100 0.208(0.036) 0.438(0.043) 0.554(0.044) n=200 0.578(0.043) 0.796(0.035) 0.870(0.029) MC15 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) 0 π 2 π 3π 2 + U2statistic n=50 0.678(0.041) 0.852(0.031) 0.908(0.025) n=100 0.894(0.027) 0.956(0.018) 0.968(0.015) n=200 0.986(0.010) 0.996(0.006) 0.998(0.004) Excess mass n=50 0.026(0.014) 0.118(0.028) 0.212(0.036) n=100 0.128(0.029) 0.318(0.041) 0.452(0.044) n=200 0.406(0.043) 0.644(0.042) 0.752(0.038) Table 3.3: Percentages of rejections for testing unimodality, with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 70 CHAPTER 3. MODE TESTING IN CIRCULAR DATA α0.01 0.05 0.10 MC11 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.006(0.007) 0.024(0.013) 0.066(0.022) n=200 0.024(0.013) 0.064(0.021) 0.090(0.025) n=1000 0.038(0.017) 0.094(0.026) 0.132(0.030) Excess mass n=50 0.010(0.009) 0.034(0.016) 0.068(0.022) n=200 0.002(0.004) 0.024(0.013) 0.058(0.020) n=1000 0(0) 0.038(0.016) 0.074(0.023) MC12 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.012(0.010) 0.048(0.019) 0.096(0.026) n=200 0.010(0.009) 0.032(0.015) 0.066(0.022) n=1000 0(0) 0.004(0.006) 0.022(0.013) Excess mass n=50 0.002(0.004) 0.028(0.014) 0.060(0.021) n=200 0.006(0.007) 0.030(0.015) 0.082(0.024) n=1000 0.004(0.006) 0.040(0.017) 0.088(0.025) MC13 0.00 0.05 0.10 0.15 0.20 0.25 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.006(0.007) 0.026(0.014) 0.058(0.020) n=200 0.002(0.004) 0.020(0.012) 0.048(0.019) n=1000 0.014(0.010) 0.044(0.018) 0.098(0.026) Excess mass n=50 0.004(0.006) 0.042(0.018) 0.100(0.026) n=200 0.008(0.008) 0.046(0.018) 0.104(0.027) n=1000 0.006(0.007) 0.056(0.020) 0.112(0.028) MC14 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.014(0.010) 0.060(0.021) 0.100(0.026) n=200 0.002(0.004) 0.020(0.012) 0.046(0.018) n=1000 0(0) 0.012(0.010) 0.030(0.015) Excess mass n=50 0.004(0.006) 0.034(0.016) 0.054(0.020) n=200 0.004(0.006) 0.036(0.016) 0.086(0.025) n=1000 0.002(0.004) 0.034(0.016) 0.082(0.024) MC15 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) 0 π 2 π 3π 2 + U2statistic n=50 0.048(0.019) 0.152(0.031) 0.238(0.037) n=200 0.040(0.017) 0.100(0.026) 0.170(0.033) n=1000 0.002(0.004) 0.036(0.016) 0.076(0.023) Excess mass n=50 0.008(0.008) 0.028(0.014) 0.066(0.022) n=200 0.012(0.010) 0.034(0.016) 0.074(0.023) n=1000 0.010(0.009) 0.040(0.017) 0.092(0.025) Table 3.4: Percentages of rejections for testing bimodality, with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 3.4. DATA ANALYSIS: DETECTION OF FIRE SEASONALITY 71 α0.01 0.05 0.10 MC16 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) 0 π 2 π 3π 2 + U2statistic n=50 0.030(0.015) 0.084(0.024) 0.156(0.032) n=200 0.026(0.014) 0.068(0.022) 0.118(0.028) n=1000 0.002(0.004) 0.010(0.009) 0.034(0.016) Excess mass n=50 0.006(0.007) 0.024(0.013) 0.060(0.021) n=200 0.002(0.004) 0.040(0.017) 0.078(0.024) n=1000 0.010(0.009) 0.048(0.019) 0.100(0.026) MC17 0.0 0.1 0.2 0.3 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.004(0.006) 0.012(0.010) 0.044(0.018) n=200 0(0) 0.012(0.010) 0.030(0.015) n=1000 0.002(0.004) 0.022(0.013) 0.052(0.019) Excess mass n=50 0.004(0.006) 0.030(0.015) 0.080(0.024) n=200 0.004(0.006) 0.036(0.016) 0.066(0.022) n=1000 0.002(0.004) 0.032(0.015) 0.072(0.023) MC18 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.002(0.004) 0.036(0.016) 0.088(0.025) n=200 0.008(0.008) 0.046(0.018) 0.126(0.029) n=1000 0.012(0.010) 0.070(0.022) 0.122(0.029) Excess mass n=50 0.010(0.009) 0.036(0.016) 0.072(0.023) n=200 0.008(0.008) 0.044(0.018) 0.074(0.023) n=1000 0.014(0.010) 0.042(0.017) 0.090(0.025) MC19 0.00 0.05 0.10 0.15 0.20 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.010(0.009) 0.044(0.018) 0.098(0.026) n=200 0.016(0.011) 0.052(0.019) 0.100(0.026) n=1000 0.016(0.011) 0.068(0.022) 0.122(0.029) Excess mass n=50 0.006(0.007) 0.046(0.018) 0.096(0.026) n=200 0.008(0.008) 0.062(0.021) 0.118(0.028) n=1000 0.012(0.010) 0.058(0.020) 0.110(0.027) MC20 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.032(0.015) 0.144(0.031) 0.214(0.036) n=200 0.054(0.020) 0.144(0.031) 0.230(0.037) n=1000 0.024(0.013) 0.102(0.027) 0.192(0.035) Excess mass n=50 0.018(0.012) 0.078(0.024) 0.140(0.030) n=200 0.016(0.011) 0.064(0.021) 0.126(0.029) n=1000 0.006(0.007) 0.040(0.017) 0.092(0.025) Table 3.5: Percentages of rejections for testing bimodality, with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 72 CHAPTER 3. MODE TESTING IN CIRCULAR DATA α0.01 0.05 0.10 MC21 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.870(0.029) 0.964(0.016) 0.978(0.013) n=100 0.998(0.004) 1(0) 1(0) n=200 1(0) 1(0) 1(0) Excess mass n=50 0.600(0.043) 0.782(0.036) 0.842(0.032) n=100 0.942(0.020) 0.978(0.013) 0.992(0.008) n=200 1(0) 1(0) 1(0) MC22 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.050(0.019) 0.180(0.034) 0.288(0.040) n=100 0.130(0.029) 0.392(0.043) 0.556(0.044) n=200 0.344(0.042) 0.700(0.040) 0.852(0.031) Excess mass n=50 0.062(0.021) 0.204(0.035) 0.282(0.039) n=100 0.184(0.034) 0.380(0.043) 0.506(0.044) n=200 0.436(0.043) 0.692(0.040) 0.780(0.036) MC23 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.032(0.015) 0.102(0.027) 0.164(0.032) n=100 0.058(0.020) 0.160(0.032) 0.218(0.036) n=200 0.060(0.021) 0.188(0.034) 0.262(0.039) Excess mass n=50 0.008(0.008) 0.054(0.020) 0.102(0.027) n=100 0.018(0.012) 0.066(0.022) 0.124(0.029) n=200 0.036(0.016) 0.088(0.025) 0.178(0.034) MC24 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π 0 π 2 π 3π 2 + U2statistic n=50 0.880(0.028) 0.962(0.017) 0.978(0.013) n=100 1(0) 1(0) 1(0) n=200 1(0) 1(0) 1(0) Excess mass n=50 0.002(0.004) 0.058(0.020) 0.108(0.027) n=100 0.022(0.013) 0.078(0.024) 0.154(0.032) n=200 0.036(0.016) 0.126(0.029) 0.200(0.035) MC25 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) 0 π 2 π 3π 2 + U2statistic n=50 0.782(0.036) 0.918(0.024) 0.950(0.019) n=100 0.980(0.012) 1(0) 1(0) n=200 0.998(0.004) 1(0) 1(0) Excess mass n=50 0.030(0.015) 0.160(0.032) 0.244(0.038) n=100 0.154(0.032) 0.340(0.042) 0.448(0.044) n=200 0.456(0.044) 0.660(0.042) 0.746(0.038) Table 3.6: Percentages of rejections for testing bimodality, with 500 simulations (1.96 times their estimated standard deviation in parenthesis) and B= 500 bootstrap samples. 3.4. DATA ANALYSIS: DETECTION OF FIRE SEASONALITY 73 3.4.1 Fires dataset The studied dataset contains the location and date of all active fires detected by the MODIS from 10 July 2002 to 9 July 2012. The MODIS algorithm identifies pixels where one or more fires are actively burning at the time of satellite overpass, based on the contrasting responses of the middle-infrared and longwave infrared bands in areas containing hot targets. Cloud and water pixels are previously excluded from analysis using multiple numerical thresholds on visible and near–infrared reflectance, and thermal infrared temperature values. The size of the smallest flaming fire having at least a 50% chance of being detected by the MODIS algorithm, under both ideal daytime and nighttime conditions is approximately 100 m2. For further details about the MODIS active fire detection algorithm, see Giglio et al. (2003) or Oom and Pereira (2013). MODIS data are provided in a discretized scale, so in order to apply the testing procedure, it is necessary to recover the continuous underlying structure. For that purpose, denote by (X1, . . . , Xn) the days of the year when the nrecorded fires occurred, with Xi∈ {1, . . . , 366}. The dataset used for the study of the number of fire seasons is the following Θi= 2π(Xi+Ei)/366; with i= 1, . . . n, being Eigenerated from the uniform distribution U(−1,0). This means that it is assumed that fires were produced at any time of the day. Provided that data are not repeated, as this can considerably alter the test statistic, other ways of modifying the data can be considered, but, in general, this perturbation does not show relevant impacts in the results. Once this modification is done, the analysed area is divided in grid cells of size 0.5◦. Then, from the resulting (360·720) cells, those ones with low fire incidence, i.e., cells having fewer than ten fires in more than 7 out of ten years, will not be considered in the study. This leaves 23,987 grid cells in the area, each one having between 37 and 20,120 fires. A map including the world map with a summary of fire counts is included in Figure 3.2. 3.4.2 Spatial False Discovered Rate As it was mentioned, removing the low fire incidence cells, the fires sample {Θi;i= 1, . . . , 4.17 ·107}is organized in 2.4·104cells. Hence, the proposed testing procedure (see Section 3.2) will be applied systematically over these groups. It is clear then 80 CHAPTER 3. MODE TESTING IN CIRCULAR DATA the ordered p–values ˘p(1) ≤. . . ≤˘p(J), unimodality is rejected in the kpatches with the smallest p–values, being k= max{υ: ˘p(υ)≤(Pυ j=1 L(j)/PJ j=1 L(j))αc} and L(j)the number of cells in the fire patch associated with ˘p(j). Step 3b. Trimming procedure. Once a decision about which patches are candidates for rejecting the null hypothesis (and hence exhibiting a multimodal fire pattern), specific cells where this rejection holds are identified. It should be noted that the cell test statistic is correlated with the test statistic at the patch level. This means that a FDR correction cannot be directly applied over all the cells belonging to the same patch and a correction is proposed by Benjamini and Hochberg (1997). First, calculate the conditional p–value of a cell within a patch that was rejected, ˆplj. Then, over these p–values, apply the two–stage procedure introduced by Benjamini et al. (2006), at level αr, to enhance the power. This last method in its first stage consists in estimating the sum of weights of null cells, using for that purpose the classical FDR procedure at level αrand then using this quantity in a second stage to determine the number of rejected cells within the patch. To be more precise and summarize Step 3b, the following steps are detailed: 4. Calculate the conditional p–value of each cell lwithin the patch that was rejected j(ˆplj). P–values can be estimated as follows: ˆplj =Z∞ zlj ˆ J0 J˜ Φ ˜ Φ−1(u1)−ˆρlju p1−ˆρ2 lj !+ + 1−ˆ J0 J!˜ Φ ˜ Φ−1(u1)−ˆρlju−ˆµj p1−ˆρ2 lj !!φ(u)du· · ˆ J0 Ju1+ 1−ˆ J0 J!˜ Φ˜ Φ−1(u1)−ˆµj!−1 , being φa standard normal density and noting that (a) u1= (Pk j=1 L(j)/PJ j=1 L(j))αcis the cutoff point of the largest p–value rejected in Step 3a. (b) ˆ J0= (J−k)/(1 −αc) is the estimated sum of weights of null patches. (c) ˆµj= ((PJ j=1 PLj l=1 zlj)/(PJ j=1 Lj))/ˆσ¯ Zjis the estimation of the standardized expectation of the patch test statistic under the alternative. 3.4. DATA ANALYSIS: DETECTION OF FIRE SEASONALITY 81 (d) ˆρlj = (1+PLj m=1,m6=lˆρj l,m)ˆσj/(Ljˆσ¯ Zj) is the estimated correlation between the z-score in a given cell and the average z-score of the patch. 5. Given these Ljp–values in the patch j, apply a two–stage procedure at level αr: (a) Apply the classic FDR procedure, at level α0 r=αr/(1 + αr). Given the ordered p–values ˆp(1)j≤. . . ≤ˆp(Lj)j, let k1j= max{υ: ˆp(υ)j≤(υ/Lj)α0 r}. (b) Apply again the classic FDR procedure at level α0 r, being in this case the sum of weights of null cells: ˆ J0j=Lj−k1j. Reject the unimodality in the k2jcells with the smallest p–values, being k2j= max{υ: ˆp(υ)j≤(υ/ ˆ J0j)α0 r}. 3.4.3 Results Employing our new proposal for testing the number of modes and correcting the FDR accounting for the spatial dependence of the data, bellow is summarized how it is applied to the wildfires dataset described in Section 3.4.1. As a first step, the p– values applying the new procedure provided in Section 3.2 (with B= 5000 bootstrap replicates) were computed in all the cells of the world. In a second step, the different fires patches were created using the land cover database. Finally, the hierarchical testing procedure was applied. First, to determine in which of the previously created patches the null hypothesis is rejected at significance level αc= 0.01. Second, within the rejected patches, it was determined which cells can be rejected at the trimming significance level αr= 0.01. The rejected cells are shown in green colour in Figure 3.5. The first main conclusion, when the modality map is observed, is that in almost all the globe, where the wildfires are present, the pattern is clearly multimodal, suggesting the prevalence of fire use as a land management tool at global scale. Although a second step indicating which is the “relevance” of each mode is needed, this map leads to the important conclusion that human are altering the fire seasonality all around the world. Although most of the non–rejected regions are scattered, one can see some areas with high–concentrated non–rejected cells. Looking at these regions, one should be careful about the conclusions related with the agricultural activity in the areas where the humans are not altering the fire seasonality. That is the case of the west coast of India, one of the largest areas with non–rejected cells, where a high number of millet crops can be found (see Leff et al., 2004), or some parts in the Central and Northern Hemisphere South America and the Indonesia islands (such us Java), as they are large producers of coffee. 82 CHAPTER 3. MODE TESTING IN CIRCULAR DATA Figure 3.5: Results after applying the procedure described in Section 3.4.2, with αc= 0.01 and αr= 0.01, in the world, divided in grids of size 0.5◦. In green: cells where H0is rejected. In blue: cells where there is no evidence to reject H0. In white, cells with low fire incidence. Chapter 4 Symmetry testing in circular data Once a conclusion about the number of modes is obtained, symmetry, or more precisely reflective symmetry, is another of these most frequently encountered simplifying assumptions. When dealing with real data it is important to see if the distribution of the sample is equal at the left and the right handsides of some central direction as it may help for a better understanding. In addition, when one is trying to figure out the underlying model of the data, the rejection of this hypothesis often leads to the subsequent exploration of models with more parameters than their symmetric counterparts. As mentioned in the Introduction, in the statistical literature, several proposals were made for testing symmetry in the linear setting. Unlike in the modality testing field, quite a few of them present a correct calibration and the issue then is trying to get the most powerful test. Related with this issue, as mentioned in the Introduction, in the linear case, some optimal, in the maximin sense (see Section 1.2.3), test have been proposed, for the two available scenarios if one wants to test symmetry: when the central location is known and when it is completely unknown. For circular data, the situation is quite different as there are just a few proposals for testing circular symmetry. Referring to the most general case, testing reflective symmetry when the central direction is unknown, just Pewsey (2002) has proposed a simple omnibus test based on the sample second sine moment about the mean direction. Following the Le Cam methodology (see Appendix D for a brief introduction to this topic), the aim of this chapter is to present new optimal tests of the null hypothesis that the circular random variable is reflectively symmetric about an unknown central direction, against the alternative hypothesis that the distribution is an asymmetric case of the k–sine–skewed family of distributions. It will be also seen that the proposal of Pewsey (2002) is (asymptotically) a particular case of the tests presented in this chapter. 83 84 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA The chapter is organized as follows. The uniform local asymptotic normality property for the k–sine–skewed family is established in Section 4.1. In Section 4.2 optimal parametric tests for circular reflective symmetry about an unknown central direction are developed, first, in Section 4.2.1, assuming that the form of the base symmetric unimodal circular density is known. Then, relaxing this assumption in Section 4.2.2 similar results are derived. In this section, the optimal semi–parametric tests for circular reflective symmetry about an unknown central direction are developed, assuming that the form of the base symmetric unimodal circular density is unknown. In Section 4.3, the results from a simulation study are presented. The simulation study is designed to explore and compare the calibration and power of the tests proposed in this chapter with those of Ley and Verdebout (2014) (when the central direction is known). In Section 4.4, the tests developed in this chapter are applied in the analysis of a circular dataset from an study about the cracks in cemented femoral components. Contents 4.1 ULAN property of k–sine–skewed densities . . . . . . . . 84 4.2 Optimal tests for reflective symmetry . . . . . . . . . . . 86 4.2.1 Optimal tests: the specified density case . . . . . . . . . . . 87 4.2.2 Optimal tests: the unspecified density case . . . . . . . . . 91 4.3 Simulation study . . . . . . . . . . . . . . . . . . . . . . . . 95 4.4 Cracks in cemented femoral components . . . . . . . . . 99 4.1 ULAN property of k–sine–skewed densities Let Θ1,...,Θnbe independent and identically distributed circular observations with common density the k–sine–skewed model (1.2.9), that is, fk µ,λ(θ) = f0(θ−µ)[1 + λsin(k(θ−µ))], where µ∈[0,2π) is the location parameter, f0∈ F a given reflectively symmetric unimodal base density and k∈N0and λ∈[−1,1] parameters controlling the number of modes and the symmetry. For any f0∈ F and any k∈N0, denote the joint distribution of the n–tuple Θ1,...,Θnby P(n) ϑ ϑ ϑ;f0,k, where ϑ ϑ ϑ= (µ, λ)0∈[0,2π)×[−1,1]. Since fk µ,λ =f0when λ= 0, and hence does not depend on k, the index kcan be dropped and simply write P(n) ϑ ϑ ϑ;f0at ϑ ϑ ϑ=ϑ ϑ ϑ0= (µ, 0)0. Any pair (f0,k) induces the parametric model P(n) f0,k=nP(n) ϑ ϑ ϑ;f0,k:ϑ ϑ ϑ∈[0,2π)×[−1,1]o, 4.1. ULAN PROPERTY OF K–SINE–SKEWED DENSITIES 85 whereas any k∈N0induces the semi–parametric model P(n) k=∪f0∈F P(n) f0,k. The uniform local asymptotic normality (ULAN) property (see Appendix D) of the parametric model P(n) f0,kin the vicinity of unimodal reflective symmetry, i.e. around λ= 0, was established by Ley and Verdebout (2014). This property of k–sine–skewed distributions is crucial in the development of the tests proposed. The derivation of the ULAN property requires the following regularity condition on the base density f0 to hold. Assumption A:The base density f0(θ)is a.e. a class–one function over [0,2π), or equivalently over Rby periodicity, with a.e.–derivative f0 0. Most classical reflectively symmetric unimodal densities, such us, the von Mises, the cardioid, the wrapped normal or the wrapped Cauchy, satisfy this requirement. Note that the continuously differentiable condition over the circle, combined with the fact that f0>0, implies that the Fisher information quantity for location, If0=R2π 0ϕ2 f0(θ)f0(θ)dθ, where ϕf0=−f0 0/f0, is finite. The ULAN property of the parametric model P(n) f0,kwith respect to ϑ ϑ ϑ= (µ, λ)0, in the vicinity of unimodal reflective symmetry, then takes the form showed in the following theorem. Theorem 4.1.1. Let f0∈ F,k∈N0and Assumption A holds. Then, for any µ∈[0,2π), the parametric family of densities P(n) f0,kis ULAN at ϑ ϑ ϑ0= (µ, 0)0with central sequence ∆ ∆ ∆(n) f0,k(µ) = ∆(n) f0,k;1(µ) ∆(n) k;2 (µ)!=1 √n n X i=1 ϕf0(Θi−µ) sin(k(Θi−µ)) , and corresponding Fisher information matrix Γ Γ Γf0,k=Γf0,k;11 Γf0,k;12 Γf0,k;12 Γf0,k;22 , where Γf0,k;11 =If0,Γf0,k;12 =−R2π 0sin(kθ)f0 0(θ)dθ and Γf0,k;22 =R2π 0sin2(kθ)f0(θ)dθ. More precisely, for any µ(n)=µ+O(n−1/2)and for any bounded sequence τ τ τ(n)= (τ(n) 1, τ(n) 2)0∈R2such that µ(n)+n−1/2τ(n) 1remains in [0,2π), then Λ(n) (µ(n)+n−1/2τ(n) 1,n−1/2τ(n) 2)0/(µ(n),0)0;f0,k= log dP(n) (µ(n)+n−1/2τ(n) 1,n−1/2τ(n) 2)0;f0,k/dP(n) (µ(n),0)0;f0 =τ τ τ(n)0∆ ∆ ∆(n) f0,k(µ(n))−(1/2)τ τ τ(n)0Γ Γ Γf0,kτ τ τ(n)+oP(1) (4.1.1) and ∆ ∆ ∆(n) f0,k(µ(n))L → N2(0 0 0,Γ Γ Γf0,k), both under P(n) (µ(n),0)0;f0as n→ ∞. 86 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA Theorem 4.1.1 is the key for providing the test statistic that will be employed and its asymptotic distribution under the null. Its proof is given in Ley and Verdebout (2014), where a brief discussion about the minimal conditions required for the ULAN property to hold are also presented. The Fisher information for departures from unimodal reflective symmetry, Γf0,k;22, and hence the cross–information quantity, Γf0,k;12, can easily be shown to be finite by bounding sin2by 1 under the integral sign. Note that the constant khas no effect on the validity of Theorem 4.1.1. For ULAN to make sense, the Fisher information matrix Γ Γ Γf0,kmust be non– singular. This is always the case, as shown in Theorem 6.1 of Ley and Verdebout (2014), except when f0is the von Mises circular density and k= 1. As it will be seen in the sequel, this singularity issue plays no role for tests for reflective symmetry about a known location but precludes the construction of a powerful test for reflective symmetry against von–Mises–based sine–skewed alternatives when the location is unknown. 4.2 Optimal tests for reflective symmetry Locally (around ϑ ϑ ϑ0) and asymptotically optimal (maximin) tests, in the Le Cam sense (see Appendix D), for reflective symmetry within the class of sine–skewed distributions which assume that µis known were proposed in Ley and Verdebout (2014). In this section, the parametric testing problem is first considered H0;f0:λ= 0 vs. Ha;f0,k:λ6= 0, for a given f0∈ F and k∈N0,(4.2.2) that is, it is studied if the underlying density of the data is f0against the alternative that it follows a k–sine–skewed version of this density, where f0is a specified density belonging to Fand the unknown central direction under H0;f0is estimated. A drawback of the above tests is that they are only valid under the parametric null hypothesis H0;f0with f0specified. In order to tackle the more general null hypothesis of reflective symmetry, a test statistic whose asymptotic distribution is valid under any symmetric density g0∈ F is needed. Thus, next the more demanding testing problem is considered H0:λ= 0 vs. Ha;k:λ6= 0, for any g0∈ F and k∈N0(4.2.3) in which both, the location parameter µand the density g0, take on nuisance roles. For both problems, the ULAN property of Theorem 4.1.1 is used to derive tests that are valid under the null hypotheses considered and achieve local and asymptotic 4.2. OPTIMAL TESTS FOR REFLECTIVE SYMMETRY 87 Location Data density Central seq. Statistic Test Section Known Known ∆(n) k;2 (µ)Q(n);µ;g0 k=|∆(n) k;2 (µ)|/Γ1/2 f0,k;22 φ(n);µ;g0 k3 of LV Unknown ∆(n) k;2 (µ)Q∗(n);µ k=|∆(n) k;2 (µ)| (Pn i=1 sin2(k(Θi−µ)))1/2φ∗(n);µ k3 of LV Unknown Known ∆(n)eff f0,k;2 (µ)Q(n);g0 kφ(n);g0 k4.2.1 Unknown ∆(n)ecd f0,g0,k;2(µ)Q∗(n) f0,kφ∗(n) f0,k4.2.2 Table 4.1: Summary of the different tests for assessing asymmetry. LV refers to Ley and Verdebout (2014); location (µ) and the data underlying density (g0∈ F) must be previously known; k∈N0and f0∈ F are selected by the user, where it is recommended to use k= 2 and VMκas f0when this information is unknown. parametric optimality against a k–sine–skewed alternative characterized by the fixed couple (f0,k)∈(F × N0). In the semi–parametric testing problem (4.2.3), f0and kare chosen a priori by the practitioner and, from these choices, tests that are asymptotically optimal against the (f0,k)–sine–skewed alternative are derived. Such tests are valid under any density g0∈ F. In order to clarify the notation in the next sections, it should be noted that ∆ denotes the central limit sequences. From them, the test statistics (denoted by Q) will be constructed, and the formal test with the asymptotic distribution (which is a standard normal for all cases) is denoted by φ. As also the Ley and Verdebout (2014) proposals will be compared in the simulation study (see Section 4.3), this will lead to a total of four different tests, reflected in Table 4.1. For Qand φ, the superindex indicates if µ∈[0,2π) or the base density of the data, g0∈ F, must be known and, with a subindex, the imposed value of k∈N0and/or the base density, not depending on the data, f0∈ F. 4.2.1 Optimal tests: the specified density case For the testing problem (4.2.2), the tests are constructed using a √nconsistent and discretized (see Assumption B below) estimator ˆµ(n). The main reason why this testing problem is more demanding than the fixed–µproblem considered in Ley and Verdebout (2014) is due to the fact that the Fisher information matrix Γ Γ Γf0,kis not, in general, diagonal. If the information matrix Γ Γ Γf0,kis diagonal, the substitution of ˆµ(n) for µwould, asymptotically, have no influence on the behavior of the central sequence 88 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA for departures from unimodal reflective symmetry ∆(n) k;2 (µ). Note that the information matrix Γ Γ Γf0,kis never diagonal if k= 1. This can be seen by noting that sin(θ)ϕf0(θ)f0(θ)>0 over (0,2π). On the other hand, when k>1 densities for which Γf0,k;12 = 0 can be found. If the density function is square integrable on [0,2π), this happens when the trigonometric moment αk=E(cos kΘ) = 0 for Θ ∼f0, which can be proved using the Fourier expansion (see Jammalamadaka and Sengupta, 2001, Section 2.1). Given any value of kwith k>1, an example where this occurs is the cardioid density for which α1=ρand αp= 0 if p > 1. According to that, the covariance Γf0,k;12 is only rarely null. Hence, a local perturbation of µhas the same asymptotic impact on ∆(n) k;2 (µ) as a local perturbation of λ= 0. It follows that the effect of ignoring the value of µis strictly positive when performing inference on λ: the stronger the correlation between µand λ, the larger that cost. The worst case occurs when the information matrix is singular (that is when f0is the von Mises and k= 1, see Section 4.1), which leads to asymptotic local powers equal to the nominal level α. This scenario implies that the best possible test is the trivial one, that is, the test that discards the observations and rejects the null hypothesis with probability α. The cost of estimating µcan be avoided by removing the effect of the location central sequence ∆(n) f0,k;1(µ) from the skewness central sequence ∆(n) k;2 (µ). To this end a Gram–Schmidt orthogonalization approach is used. This is done projecting ∆(n) k;2 (µ) onto the subspace orthogonal to ∆(n) f0,k;1(µ), which ensures that the resulting f0– efficient central sequence for skewness ∆(n)eff f0,k;2 (µ) and ∆(n) f0,k;1(µ) are asymptotically uncorrelated. This new central sequence is thus of the form ∆(n)eff f0,k;2 (µ) = ∆(n) k;2 (µ)−Γf0,k;12 Γf0,k;11 ∆(n) f0,k;1(µ) =n−1/2 n X i=1 sin(k(Θi−µ)) −Γf0,k;12 Γf0,k;11 ϕf0(Θi−µ).(4.2.4) At this point, another important consequence of the ULAN property should be employed, the asymptotic linearity property (see Appendix D): ∆(n) f0,k(µ+n−1/2τ(n) 1)−∆(n) f0,k(µ) = −Γf0,k(τ(n) 1,0)0+oP(1) (4.2.5) under P(n) (µ,0)0;f0as n→ ∞, with τ(n) 1∈Ras in Theorem 4.1.1. See Sections 2 and 3 of Koudou and Ley (2014) for further explanations on these issues. It is now not difficult to derive from (4.2.5) the asymptotic linearity property of ∆(n)eff f0,k;2 (µ): ∆(n)eff f0,k;2 (µ+n−1/2τ(n) 1)−∆(n)eff f0,k;2 (µ) = oP(1) (4.2.6) 4.2. OPTIMAL TESTS FOR REFLECTIVE SYMMETRY 89 under P(n) (µ,0)0;f0as n→ ∞. Now consider replacing the non–random bounded sequence τ(n) 1with n1/2(ˆµ(n)−µ) for some √nconsistent estimator ˆµ(n). The latter is bounded in probability and, via Lemma 4.4 of Kreiss (1987), serves as a candidate for τ(n) 1 provided the following assumption holds. Assumption B:The sequence of estimators ˆµ(n)is (i) √nconsistent, i.e., n1/2(ˆµ(n) −µ) = OP(1) as n→ ∞, under P(n) (µ,0)0;f0, and (ii) locally asymptotically discrete, meaning that, for all µ∈[0,2π)and all c > 0, there exists an M=M(c)>0such that the number of possible values of ˆµ(n)in intervals of the form {t∈R:n1/2|t−µ| ≤ c} is bounded by M, uniformly as n→ ∞. Note that Assumption B (ii) is a purely technical requirement, with little practical implication. Indeed, for fixed sample size, any estimator can be considered part of a locally asymptotically discrete sequence. However, it is precisely this assumption that allows us to replace τ(n) 1by n1/2(ˆµ(n)−µ) in (4.2.6) thanks to the aforementioned Lemma 4.4 of Kreiss (1987), yielding ∆(n)eff f0,k;2 (ˆµ(n))−∆(n)eff f0,k;2 (µ) = oP(1) (4.2.7) under P(n) (µ,0)0;f0as n→ ∞. The locally (around λ= 0) and asymptotically maximin f0–parametric test φ(n);f0 kcan be built. It rejects H0;f0at asymptotic level αwhenever the statistic Q(n);f0 k=|∆(n)eff f0,k;2 (ˆµ(n))| Γ1/2 f0,k;22.1 exceeds the upper α/2 quantile of the standard normal distribution, z1−α/2, where Γf0,k;22.1= Γf0,k;22 −Γ2 f0,k;12/Γf0,k;11 is the asymptotic variance of ∆(n)eff f0,k;2 (µ) under P(n) (µ,0)0;f0. Optimal properties of this test statistic will be described in a more flexible scenario in Section 4.2.2. Depending on the choice of f0and k, there are different constructions of the test statistic Q(n);f0 k. Among the possible base symmetric densities, here the test statistic is described for three of the most well–known distributions: the von Mises, the cardioid and the wrapped Cauchy. Von Mises distribution For the von Mises distribution, it is obtained that ϕfVMκ(θ) = κsin(θ), ΓfVMκ,k;11 = κI1(κ)/I0(κ), ΓfVMκ,k;12 =kIk(κ)/I0(κ) and ΓfVMκ,k;22 = (1 −I2k(κ)/I0(κ))/2, where Ikdenotes the modified Bessel function of the first kind and order k. As it was 96 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) π3π20π2π 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) π3π20π2π 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) π3π20π2π 0.0 0.2 0.4 0.6 0.8 1.0 1.2 θ f(θ) π3π20π2π 0.0 0.2 0.4 0.6 0.8 1.0 1.2 θ f(θ) π3π20π2π 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 θ f(θ) π3π20π2π Figure 4.1: Examples of k–sine–skewed distributions. Left: λ= 0 (symmetric). Center: λ= 0.2 and k= 1 (“low” asymmetry). Right: λ= 0.6 and k= 3 (“high” asymmetry). First row: fVM1(solid line), fWN0.5(dashed line) and fWC0.5(dotted line). Second row: fVM10 (solid line), fWN0.9(dashed line) and fC0.45 (dotted line). φ(n);fVMκ 1, the resulting test coincides with the trivial test and for that reason is not showed on the simulation results in Table 4.5. General conclusions. First, as mentioned in the previous section, two decisions must be done, selecting the base density f0(in the semiparametric test for unknown location) and the value of k. Broadly, as it will be seen later in this simulation study, when the value of kis unknown the best power results are achieved when k= 2. Although it will be seen that, in practice, the selection of the base density f0in the semiparametric scenario has little impact on the test statistic, our suggestion is employ the von Mises density as, independently on the concentration parameter, it is optimal against all the kind of von-Mises and cardioid 2–sine–skewed alternatives. From the results reported in Tables 4.2–4.4 it may seem that the calibration of the tests is correct; although for the scenarios in which the value of µis unspecified, the φ(n);g0 2and φ∗(n) fVMκ;2 tests are somewhat conservative for n= 30. As might be expected, the rejection rates of all four tests increase with the sample size, n, and 4.3. SIMULATION STUDY 97 the value of λ, and are generally higher when the posited value of kcoincides with the true value of k0. In general, knowing the form of the underlying distribution has little effect on the relative performance of the two pairs of tests. Again as might be expected, the rejection rates are higher for the tests for which the centre, µ, is correctly imposed. For those two test, φ(n);g0 kand φ∗(n) fVMκ;k, the power is lowest when the underlying distribution is highly concentrated, as it is for the fVM10 and fWN0.9 cases. The reason for this lower power can be appreciated from a comparison of their densities in the graphics showed in Figure 4.1. Clearly, there is very little difference between the shape of the base symmetric density and the shapes of the other densities, corresponding to different values of λ, portrayed in those graphics. Overall, however, φ∗(n) fVMκ;kprovides a relatively competitive test compared with φ(n);g0 kfor the scenario in which µis unspecified. Selection of k. In this part of the simulation study, the four tests for reflective symmetry assuming that kwas 1, 2 and 3 and imposing the cardioid C0.25 as the base density for the semi–parametric test, φ∗(n) fC0.25 ;k. This base density was chosen as, when k= 2, it presents similar results to those achieved with the von Mises base density and, unlike this density, it also allows to impose k= 1 in the construction of the test. The rejection rates obtained for the different values of kand a nominal significance level of α= 0.05 are reproduced in Tables 4.5 and 4.6 (for k= 1), Tables 4.2–4.4 (for k= 2) and in Tables 4.7 and 4.8 (for k= 3). From the results reported in Tables 4.2–4.8 it may seem that the rejection rates of all four tests are generally highest when the posited value of kcoincides with the true value of k0. Some exceptions to this rule can be found in the concentrated densities fVM10 and fWN0.9. For example, for these two densities, in Tables 4.2–4.6, when k= 1,2; the φ(n);µ;g0 kand φ∗(n);µ ktests perform better for k0= 3, rather than for k0=k. When k6=k0, some of the tests perform, at best, like the trivial test. This is the case, in Tables 4.2–4.8, for the φ(n);g0 kand φ∗(n) fC0.25 ;ktests when k0= 1 and, again, g0 is the highly concentrated fVM10 or fWN0.9density. This last case also occurs for other densities as, for example, the φ(n);µ;g0 1and φ∗(n);µ 1tests when g0=fC0.45 and k0= 3, in Table 4.5; φ(n);µ;g0 3,φ∗(n);µ 3,φ(n);g0 3and φ∗(n) fC0.25 ;3 tests when g0=fC0.45 and k0= 1, in Table 4.7; φ(n);g0 3and φ∗(n) fC0.25 ;3 tests when g0=fWN0.5or g0=fWN0.9and k0= 1, in Table 4.8. Using the asymptotic distribution for calibrating the test. To improve the calibration results when the sample size is too small (n= 30), one can think in nonparametric techniques. Replicating the proposal of Pewsey (2002), the employed method, for generating resamples under the null hypothesis of symmetry, was sampling with replacement from {Θ1,...,Θn,2ˆµ(n)−Θ1,2ˆµ(n)−Θn}to generate the 98 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA new bootstrap resamples of size nand, then, compute their associated test statistic value. Observing the obtained results in Table 4.9, the conclusion is that it is not worthwhile the extra computational time to obtain the p–value with the bootstrap method. Wrong selection of the base density in the parametric test. In order to determine the behaviour of the parametric tests for reflective symmetry under different densities from those used as base symmetric densities, in Tables 4.10, 4.11 and 4.12, the parametric test (with k= 2) was used for computing the percentage of rejections, using as f0the following densities: fVM10 ,fCρand fWC0.5. In Tables 4.10– 4.12, it is observed the high impact of wrongly determine the base density on the parametric test statistic, specially when the location parameter is unknown, a clear example of this fact can be found (in Table 4.12) in the poor calibration of the parametric test (using the wrong underlying density) when the sample is generated from the wrapped normal density with concentration parameter 0.9. Using other base densities in the semiparametric test. Tables 4.13 and 4.14 together with Tables 4.2–4.4 can be used to compare the percentage of rejections of the semiparametric proposal with unknown central direction (and k= 2) under other base densities f0:fC0.45 and fWC0.5. From Tables 4.13 and 4.14, the little impact of the base density in the semiparametric scenario is observed, since they all provide, in the reported scenarios, a similar percentage of rejections as in Tables 4.2–4.4. Using the Ley and Verdebout (2014) proposal when µis unknown. Also, Tables 4.13 and 4.14 report the percentage of rejections when it is used the semiparametric proposal with “known” central direction estimated from the sample, φ∗(n);ˆµ(n) 2, in order to analyse the impact of estimating µin the Ley and Verdebout (2014) proposal. On these tables, the already mentioned high impact of using the test with known location parameter, when it is unknown (estimating it), is observed. A clear example of the bad calibration using the Ley and Verdebout (2014) proposal estimating the location parameter from the sample can be found when the data is generated from concentrated densities as fWN0.9(Table 4.17). In Tables 4.15 and 4.16 the behaviour of the known location test when imposing a wrong parameter of µcan be observed. In most of the cases, even an slight perturbation in the location parameter has an important effect in the calibration of the test. For example, imposing a wrong location, translated by π/16, for large sample sizes (n= 500) in all the simulated densities, with the exception of fC0.45 , the percentage of rejections under the null is always above the significance level. Suggesting that the Ley and Verdebout (2014) proposal should not be used unless the location parameter is clearly known. 4.4. CRACKS IN CEMENTED FEMORAL COMPONENTS 99 ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ●● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ●● ● ● ● ●● ● ● ●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ●● ● ● ●● ● ●● ● ●● ● ● ● ● ●● ● ● ● ● ●● ● ● ● ● ● ●● ● ● ● ● ●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ●● ● ● ● ● ●● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● Lateral Anterior Medial Posterior + ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ●● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●● ● ● ● ● ●● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● Lateral Anterior Medial Posterior + Figure 4.2: Raw circular plot, rose diagram and density estimate for the cracks in the cemented femoral components in the proximal (left) and the distal region (right). Behaviour of the test under other models. In Table 4.17, the percentage of rejections using the semiparametric test with known and unknown central direction are calculated under other alternatives. With this objective, it was simulated samples of size n= 30,100,500 from the Kato and Jones (2012) distribution, with parameters µ= 0, r= 0.5, κ= 0.5,0.9 and skewness parameter ν= 0,0.2,0.4,0.6; and the three–parameter asymmetric submodel given in the Equation (7) of Kato and Jones (2015), with parameters µ= 0, γ= 0.5,0.9, ¯ β2=νγ(1 −γ) and skewness parameter ν= 0,0.2,0.4,0.6 (see Appendix A). Under these other alternatives, the test with known and unknown central location provides, in general, good results as the power increases with the sample size and the skewness of the model. 4.4 Cracks in cemented femoral components To illustrate the application of the new test for circular symmetry, the angles of the fatigue cracks around the cemented femoral components, described in Section 1.3.4, will be considered. The objective is to test if these cracks are symmetrically distributed around a central direction in both regions distinguished in the femur (proximal, close to the hip, and distal, close to the knee). In this study the entire sample will be taken, with the exception of an outlier bone, as Mann et al. (2003) reported that this specimen had an inferior cement mantle with substantial stem–cement voids and, 100 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 k0g0=fVM1 Q(n);µ;g0 21 0.050 0.046 0.050 0.057 0.087 0.306 0.097 0.232 0.797 0.145 0.463 0.990 Q∗(n);µ 20.053 0.043 0.048 0.051 0.088 0.298 0.095 0.232 0.802 0.148 0.473 0.990 Q(n);g0 20.057 0.043 0.056 0.046 0.046 0.054 0.046 0.044 0.068 0.036 0.064 0.219 Q∗(n) fVMκ;2 0.033 0.036 0.050 0.029 0.043 0.056 0.038 0.038 0.068 0.044 0.079 0.248 Q∗(n) fC0.25 ;2 0.033 0.037 0.048 0.030 0.045 0.056 0.039 0.040 0.066 0.045 0.078 0.250 Q(n);µ;g0 22 0.050 0.046 0.050 0.100 0.295 0.895 0.349 0.825 1 0.648 0.994 1 Q∗(n);µ 20.053 0.043 0.048 0.111 0.293 0.897 0.347 0.826 1 0.646 0.995 1 Q(n);g0 20.057 0.043 0.056 0.069 0.211 0.783 0.195 0.605 0.999 0.351 0.826 1 Q∗(n) fVMκ;2 0.033 0.036 0.050 0.039 0.183 0.777 0.139 0.550 0.999 0.235 0.747 1 Q∗(n) fC0.25 ;2 0.033 0.037 0.048 0.043 0.180 0.773 0.140 0.551 0.999 0.237 0.758 1 Q(n);µ;g0 23 0.050 0.046 0.050 0.062 0.100 0.321 0.103 0.256 0.823 0.164 0.475 0.991 Q∗(n);µ 20.053 0.043 0.048 0.061 0.105 0.319 0.107 0.259 0.823 0.172 0.482 0.991 Q(n);g0 20.057 0.043 0.056 0.053 0.077 0.293 0.086 0.218 0.794 0.122 0.405 0.982 Q∗(n) fVMκ;2 0.033 0.036 0.050 0.032 0.065 0.288 0.057 0.198 0.777 0.086 0.343 0.975 Q∗(n) fC0.25 ;2 0.033 0.037 0.048 0.029 0.059 0.262 0.058 0.182 0.725 0.082 0.306 0.952 k0g0=fVM10 Q(n);µ;g0 21 0.049 0.046 0.044 0.060 0.085 0.282 0.097 0.233 0.786 0.163 0.476 0.983 Q∗(n);µ 20.051 0.042 0.044 0.058 0.090 0.282 0.095 0.234 0.791 0.162 0.464 0.986 Q(n);g0 20.050 0.047 0.040 0.045 0.045 0.034 0.045 0.040 0.039 0.047 0.036 0.044 Q∗(n) fVMκ;2 0.028 0.040 0.040 0.029 0.037 0.034 0.029 0.030 0.041 0.027 0.038 0.053 Q∗(n) fC0.25 ;2 0.028 0.041 0.039 0.029 0.037 0.034 0.028 0.031 0.041 0.027 0.038 0.053 Q(n);µ;g0 22 0.049 0.046 0.044 0.083 0.183 0.694 0.222 0.586 0.998 0.417 0.911 1 Q∗(n);µ 20.051 0.042 0.044 0.083 0.190 0.699 0.230 0.586 0.998 0.431 0.918 1 Q(n);g0 20.050 0.047 0.040 0.052 0.053 0.044 0.045 0.051 0.070 0.032 0.040 0.053 Q∗(n) fVMκ;2 0.028 0.040 0.040 0.029 0.044 0.046 0.038 0.042 0.068 0.037 0.047 0.068 Q∗(n) fC0.25 ;2 0.028 0.041 0.039 0.028 0.044 0.045 0.037 0.042 0.068 0.035 0.047 0.067 Q(n);µ;g0 23 0.049 0.046 0.044 0.107 0.238 0.843 0.273 0.728 1 0.575 0.977 1 Q∗(n);µ 20.051 0.042 0.044 0.104 0.243 0.842 0.280 0.743 1 0.581 0.982 1 Q(n);g0 20.050 0.047 0.040 0.061 0.058 0.138 0.050 0.098 0.412 0.059 0.138 0.543 Q∗(n) fVMκ;2 0.028 0.040 0.040 0.034 0.057 0.146 0.033 0.087 0.397 0.031 0.091 0.524 Q∗(n) fC0.25 ;2 0.028 0.041 0.039 0.033 0.057 0.144 0.033 0.086 0.400 0.030 0.089 0.524 Table 4.2: Percentages of rejections for testing symmetry, with 1000 simulations. 4.4. CRACKS IN CEMENTED FEMORAL COMPONENTS 101 λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 k0g0=fC0.45 Q(n);µ;g0 21 0.056 0.057 0.042 0.068 0.090 0.285 0.106 0.260 0.814 0.175 0.489 0.993 Q∗(n);µ 20.061 0.057 0.041 0.070 0.093 0.287 0.103 0.252 0.811 0.171 0.484 0.992 Q(n);g0 20.054 0.061 0.044 0.041 0.064 0.246 0.042 0.101 0.485 0.025 0.055 0.338 Q∗(n) fVMκ;2 0.049 0.057 0.038 0.045 0.067 0.269 0.055 0.148 0.577 0.048 0.118 0.488 Q∗(n) fC0.25 ;2 0.049 0.057 0.038 0.044 0.067 0.272 0.056 0.150 0.586 0.048 0.120 0.504 Q(n);µ;g0 22 0.056 0.057 0.042 0.120 0.289 0.905 0.336 0.830 1 0.661 0.991 1 Q∗(n);µ 20.061 0.057 0.041 0.119 0.292 0.905 0.342 0.834 1 0.667 0.991 1 Q(n);g0 20.054 0.061 0.044 0.073 0.261 0.890 0.200 0.729 1 0.402 0.948 1 Q∗(n) fVMκ;2 0.049 0.057 0.038 0.083 0.274 0.901 0.205 0.763 1 0.373 0.950 1 Q∗(n) fC0.25 ;2 0.049 0.057 0.038 0.084 0.275 0.901 0.210 0.767 1 0.378 0.957 1 Q(n);µ;g0 23 0.056 0.057 0.042 0.057 0.094 0.296 0.113 0.237 0.829 0.157 0.493 0.989 Q∗(n);µ 20.061 0.057 0.041 0.058 0.100 0.299 0.112 0.237 0.832 0.157 0.496 0.988 Q(n);g0 20.054 0.061 0.044 0.051 0.092 0.287 0.070 0.217 0.816 0.105 0.435 0.987 Q∗(n) fVMκ;2 0.049 0.057 0.038 0.050 0.085 0.302 0.074 0.215 0.802 0.106 0.408 0.981 Q∗(n) fC0.25 ;2 0.049 0.057 0.038 0.046 0.086 0.302 0.075 0.213 0.806 0.107 0.399 0.984 k0g0=fWN0.5 Q(n);µ;g0 21 0.054 0.048 0.052 0.066 0.113 0.354 0.117 0.292 0.894 0.202 0.578 0.997 Q∗(n);µ 20.052 0.053 0.049 0.065 0.112 0.351 0.111 0.294 0.895 0.212 0.579 0.998 Q(n);g0 20.036 0.047 0.052 0.040 0.059 0.145 0.033 0.066 0.223 0.021 0.040 0.112 Q∗(n) fVMκ;2 0.030 0.041 0.050 0.039 0.059 0.148 0.037 0.077 0.257 0.041 0.058 0.167 Q∗(n) fC0.25 ;2 0.028 0.042 0.049 0.040 0.060 0.149 0.036 0.077 0.268 0.043 0.062 0.170 Q(n);µ;g0 22 0.054 0.048 0.052 0.116 0.289 0.894 0.327 0.826 1 0.656 0.995 1 Q∗(n);µ 20.052 0.053 0.049 0.112 0.289 0.896 0.329 0.825 1 0.669 0.996 1 Q(n);g0 20.036 0.047 0.052 0.090 0.243 0.864 0.210 0.717 1 0.399 0.941 1 Q∗(n) fVMκ;2 0.030 0.041 0.050 0.084 0.249 0.866 0.186 0.728 1 0.355 0.930 1 Q∗(n) fC0.25 ;2 0.028 0.042 0.049 0.085 0.254 0.869 0.187 0.731 1 0.359 0.935 1 Q(n);µ;g0 23 0.054 0.048 0.052 0.062 0.118 0.364 0.114 0.284 0.889 0.187 0.553 0.999 Q∗(n);µ 20.052 0.053 0.049 0.061 0.118 0.363 0.114 0.277 0.889 0.189 0.558 0.999 Q(n);g0 20.036 0.047 0.052 0.050 0.112 0.354 0.084 0.265 0.893 0.167 0.531 0.999 Q∗(n) fVMκ;2 0.030 0.041 0.050 0.062 0.111 0.387 0.086 0.269 0.912 0.165 0.534 1 Q∗(n) fC0.25 ;2 0.028 0.042 0.049 0.062 0.108 0.374 0.082 0.259 0.903 0.152 0.509 1 Table 4.3: Percentages of rejections for testing symmetry, with 1000 simulations. 102 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 k0g0=fWN0.9 Q(n);µ;g0 21 0.043 0.047 0.053 0.059 0.109 0.437 0.109 0.339 0.954 0.227 0.671 1 Q∗(n);µ 20.050 0.041 0.055 0.069 0.112 0.434 0.100 0.340 0.953 0.230 0.679 1 Q(n);g0 20.052 0.050 0.059 0.044 0.047 0.051 0.047 0.041 0.058 0.030 0.031 0.051 Q∗(n) fVMκ;2 0.024 0.044 0.060 0.033 0.044 0.055 0.032 0.039 0.063 0.033 0.036 0.065 Q∗(n) fC0.25 ;2 0.024 0.046 0.057 0.032 0.044 0.055 0.028 0.039 0.065 0.032 0.035 0.063 Q(n);µ;g0 22 0.043 0.047 0.053 0.089 0.219 0.824 0.263 0.717 1 0.564 0.982 1 Q∗(n);µ 20.050 0.041 0.055 0.091 0.227 0.825 0.263 0.731 1 0.571 0.983 1 Q(n);g0 20.052 0.050 0.059 0.050 0.068 0.128 0.055 0.088 0.354 0.051 0.124 0.446 Q∗(n) fVMκ;2 0.024 0.044 0.060 0.043 0.051 0.133 0.031 0.066 0.357 0.024 0.104 0.462 Q∗(n) fC0.25 ;2 0.024 0.046 0.057 0.042 0.050 0.129 0.030 0.063 0.359 0.023 0.098 0.463 Q(n);µ;g0 23 0.043 0.047 0.053 0.097 0.233 0.829 0.275 0.737 1 0.599 0.986 1 Q∗(n);µ 20.050 0.041 0.055 0.096 0.232 0.824 0.280 0.741 1 0.608 0.988 1 Q(n);g0 20.052 0.050 0.059 0.073 0.113 0.483 0.131 0.326 0.951 0.194 0.546 0.998 Q∗(n) fVMκ;2 0.024 0.044 0.060 0.049 0.123 0.486 0.077 0.297 0.958 0.123 0.497 0.997 Q∗(n) fC0.25 ;2 0.024 0.046 0.057 0.048 0.120 0.483 0.076 0.291 0.956 0.122 0.493 0.997 k0g0=fWC0.5 Q(n);µ;g0 21 0.058 0.041 0.045 0.061 0.098 0.241 0.092 0.222 0.682 0.132 0.389 0.959 Q∗(n);µ 20.054 0.043 0.047 0.062 0.097 0.241 0.084 0.213 0.691 0.133 0.387 0.962 Q(n);g0 20.045 0.054 0.043 0.046 0.064 0.162 0.059 0.122 0.507 0.083 0.264 0.894 Q∗(n) fVMκ;2 0.038 0.053 0.053 0.051 0.069 0.227 0.077 0.187 0.719 0.120 0.427 0.979 Q∗(n) fC0.25 ;2 0.038 0.053 0.050 0.051 0.071 0.220 0.076 0.187 0.703 0.117 0.412 0.974 Q(n);µ;g0 22 0.058 0.041 0.045 0.114 0.294 0.852 0.322 0.807 1 0.629 0.985 1 Q∗(n);µ 20.054 0.043 0.047 0.110 0.292 0.852 0.325 0.808 1 0.624 0.987 1 Q(n);g0 20.045 0.054 0.043 0.069 0.144 0.478 0.113 0.372 0.957 0.219 0.628 0.997 Q∗(n) fVMκ;2 0.038 0.053 0.053 0.060 0.123 0.417 0.106 0.310 0.910 0.160 0.467 0.990 Q∗(n) fC0.25 ;2 0.038 0.053 0.050 0.064 0.130 0.432 0.107 0.320 0.914 0.168 0.485 0.991 Q(n);µ;g0 23 0.058 0.041 0.045 0.067 0.113 0.334 0.110 0.290 0.857 0.186 0.531 0.997 Q∗(n);µ 20.054 0.043 0.047 0.067 0.114 0.336 0.110 0.289 0.859 0.190 0.543 0.996 Q(n);g0 20.045 0.054 0.043 0.046 0.061 0.057 0.056 0.060 0.088 0.047 0.060 0.127 Q∗(n) fVMκ;2 0.038 0.053 0.053 0.051 0.074 0.152 0.074 0.128 0.413 0.097 0.206 0.679 Q∗(n) fC0.25 ;2 0.038 0.053 0.050 0.049 0.072 0.127 0.072 0.111 0.323 0.088 0.169 0.552 Table 4.4: Percentages of rejections for testing symmetry, with 1000 simulations. 4.4. CRACKS IN CEMENTED FEMORAL COMPONENTS 103 λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 k0g0=fVM1 Q(n);µ;g0 11 0.040 0.050 0.052 0.125 0.302 0.862 0.326 0.788 1 0.636 0.995 1 Q∗(n);µ 10.043 0.052 0.052 0.115 0.302 0.865 0.330 0.791 1 0.650 0.996 1 Q∗(n) fC0.25 ;1 0.033 0.036 0.045 0.035 0.043 0.052 0.045 0.040 0.065 0.046 0.073 0.245 Q(n);µ;g0 12 0.040 0.050 0.052 0.076 0.116 0.337 0.127 0.281 0.845 0.196 0.512 0.993 Q∗(n);µ 10.043 0.052 0.052 0.074 0.114 0.334 0.127 0.279 0.845 0.190 0.516 0.992 Q∗(n) fC0.25 ;1 0.033 0.036 0.045 0.045 0.177 0.740 0.133 0.535 0.998 0.220 0.753 1 Q(n);µ;g0 13 0.040 0.050 0.052 0.058 0.061 0.084 0.052 0.051 0.131 0.060 0.083 0.204 Q∗(n);µ 10.043 0.052 0.052 0.060 0.060 0.085 0.049 0.057 0.127 0.061 0.084 0.203 Q∗(n) fC0.25 ;1 0.033 0.036 0.045 0.029 0.050 0.105 0.044 0.082 0.263 0.042 0.109 0.462 k0g0=fVM10 Q(n);µ;g0 11 0.051 0.046 0.043 0.059 0.093 0.282 0.101 0.241 0.788 0.168 0.491 0.984 Q∗(n);µ 10.045 0.043 0.044 0.057 0.095 0.285 0.101 0.235 0.793 0.167 0.482 0.986 Q∗(n) fC0.25 ;1 0.027 0.039 0.038 0.028 0.036 0.035 0.028 0.031 0.041 0.027 0.037 0.052 Q(n);µ;g0 12 0.051 0.046 0.043 0.082 0.179 0.681 0.224 0.575 0.998 0.413 0.903 1 Q∗(n);µ 10.045 0.043 0.044 0.078 0.186 0.693 0.218 0.584 0.997 0.434 0.913 1 Q∗(n) fC0.25 ;1 0.027 0.039 0.038 0.028 0.043 0.043 0.034 0.041 0.068 0.034 0.046 0.065 Q(n);µ;g0 13 0.051 0.046 0.043 0.101 0.221 0.811 0.259 0.702 1 0.535 0.970 1 Q∗(n);µ 10.045 0.043 0.044 0.101 0.231 0.808 0.263 0.718 1 0.563 0.975 1 Q∗(n) fC0.25 ;1 0.027 0.039 0.038 0.032 0.058 0.142 0.033 0.078 0.393 0.030 0.088 0.521 k0g0=fC0.45 Q(n);µ;g0 11 0.039 0.053 0.057 0.123 0.299 0.876 0.334 0.829 1 0.673 0.996 1 Q∗(n);µ 10.036 0.050 0.056 0.121 0.306 0.880 0.332 0.834 1 0.679 0.998 1 Q(n);g0 10.059 0.052 0.052 0.046 0.087 0.311 0.048 0.133 0.630 0.035 0.069 0.499 Q∗(n) fC0.25 ;1 0.048 0.057 0.046 0.048 0.088 0.311 0.052 0.180 0.685 0.051 0.162 0.648 Q(n);µ;g0 12 0.039 0.053 0.057 0.063 0.088 0.311 0.105 0.239 0.812 0.179 0.461 0.988 Q∗(n);µ 10.036 0.050 0.056 0.059 0.091 0.311 0.108 0.242 0.811 0.177 0.466 0.989 Q(n);g0 10.059 0.052 0.052 0.083 0.193 0.764 0.169 0.588 1 0.321 0.903 1 Q∗(n) fC0.25 ;1 0.048 0.057 0.046 0.086 0.264 0.884 0.209 0.750 1 0.366 0.958 1 Q(n);µ;g0 13 0.039 0.053 0.057 0.057 0.048 0.064 0.067 0.042 0.053 0.059 0.052 0.038 Q∗(n);µ 10.036 0.050 0.056 0.055 0.048 0.063 0.068 0.045 0.052 0.064 0.048 0.036 Q(n);g0 10.059 0.052 0.052 0.048 0.054 0.066 0.040 0.045 0.037 0.039 0.047 0.045 Q∗(n) fC0.25 ;1 0.048 0.057 0.046 0.046 0.068 0.123 0.056 0.085 0.327 0.057 0.178 0.605 Table 4.5: Percentages of rejections for testing symmetry, with 1000 simulations. 104 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 k0g0=fWN0.5 Q(n);µ;g0 11 0.045 0.048 0.046 0.114 0.238 0.879 0.301 0.792 1 0.629 0.990 1 Q∗(n);µ 10.046 0.049 0.047 0.109 0.233 0.880 0.303 0.789 1 0.647 0.991 1 Q(n);g0 10.039 0.047 0.063 0.052 0.048 0.155 0.041 0.065 0.276 0.022 0.044 0.152 Q∗(n) fC0.25 ;1 0.019 0.048 0.058 0.043 0.056 0.168 0.033 0.090 0.312 0.036 0.068 0.222 Q(n);µ;g0 12 0.045 0.048 0.046 0.052 0.113 0.357 0.119 0.296 0.913 0.209 0.586 0.999 Q∗(n);µ 10.046 0.049 0.047 0.057 0.112 0.354 0.120 0.301 0.914 0.207 0.588 0.999 Q(n);g0 10.039 0.047 0.063 0.083 0.218 0.783 0.176 0.626 0.998 0.361 0.902 1 Q∗(n) fC0.25 ;1 0.019 0.048 0.058 0.080 0.243 0.851 0.177 0.714 1 0.353 0.939 1 Q(n);µ;g0 13 0.045 0.048 0.046 0.049 0.043 0.057 0.051 0.048 0.088 0.041 0.049 0.109 Q∗(n);µ 10.046 0.049 0.047 0.051 0.045 0.060 0.050 0.049 0.089 0.047 0.052 0.107 Q(n);g0 10.039 0.047 0.063 0.050 0.064 0.088 0.047 0.070 0.186 0.044 0.100 0.353 Q∗(n) fC0.25 ;1 0.019 0.048 0.058 0.058 0.084 0.169 0.052 0.140 0.525 0.095 0.237 0.866 k0g0=fWN0.9 Q(n);µ;g0 11 0.037 0.042 0.055 0.064 0.117 0.460 0.111 0.367 0.966 0.236 0.694 1 Q∗(n);µ 10.048 0.043 0.054 0.067 0.115 0.462 0.113 0.364 0.967 0.253 0.703 1 Q(n);g0 10.051 0.056 0.052 0.045 0.049 0.052 0.043 0.041 0.062 0.032 0.030 0.045 Q∗(n) fC0.25 ;1 0.023 0.043 0.057 0.029 0.038 0.053 0.020 0.037 0.061 0.029 0.031 0.059 Q(n);µ;g0 12 0.037 0.042 0.055 0.098 0.213 0.787 0.239 0.685 1 0.532 0.974 1 Q∗(n);µ 10.048 0.043 0.054 0.093 0.212 0.789 0.249 0.703 1 0.558 0.979 1 Q(n);g0 10.051 0.056 0.052 0.048 0.068 0.133 0.054 0.087 0.349 0.055 0.140 0.460 Q∗(n) fC0.25 ;1 0.023 0.043 0.057 0.037 0.046 0.126 0.028 0.057 0.355 0.020 0.092 0.470 Q(n);µ;g0 13 0.037 0.042 0.055 0.085 0.183 0.705 0.223 0.637 1 0.470 0.946 1 Q∗(n);µ 10.048 0.043 0.054 0.081 0.191 0.708 0.237 0.640 1 0.505 0.947 1 Q(n);g0 10.051 0.056 0.052 0.073 0.108 0.446 0.128 0.308 0.941 0.189 0.519 0.998 Q∗(n) fC0.25 ;1 0.023 0.043 0.057 0.041 0.104 0.470 0.071 0.275 0.948 0.111 0.480 0.995 k0g0=fWC0.5 Q(n);µ;g0 11 0.051 0.056 0.046 0.098 0.241 0.763 0.268 0.708 1 0.544 0.963 1 Q∗(n);µ 10.050 0.055 0.044 0.092 0.231 0.763 0.273 0.716 1 0.561 0.967 1 Q(n);g0 10.037 0.048 0.052 0.041 0.072 0.263 0.062 0.173 0.742 0.102 0.363 0.979 Q∗(n) fC0.25 ;1 0.036 0.056 0.045 0.050 0.068 0.198 0.062 0.163 0.616 0.098 0.358 0.949 Q(n);µ;g0 12 0.051 0.056 0.046 0.066 0.083 0.257 0.103 0.233 0.769 0.165 0.454 0.982 Q∗(n);µ 10.050 0.055 0.044 0.061 0.090 0.259 0.114 0.227 0.768 0.171 0.456 0.981 Q(n);g0 10.037 0.048 0.052 0.047 0.077 0.271 0.090 0.217 0.707 0.137 0.344 0.892 Q∗(n) fC0.25 ;1 0.036 0.056 0.045 0.063 0.132 0.450 0.100 0.324 0.923 0.163 0.506 0.994 Q(n);µ;g0 13 0.051 0.056 0.046 0.055 0.057 0.093 0.046 0.110 0.256 0.067 0.148 0.534 Q∗(n);µ 10.050 0.055 0.044 0.055 0.051 0.097 0.054 0.110 0.262 0.076 0.152 0.533 Q(n);g0 10.037 0.048 0.052 0.047 0.119 0.392 0.111 0.326 0.897 0.184 0.546 0.990 Q∗(n) fC0.25 ;1 0.036 0.056 0.045 0.041 0.064 0.057 0.056 0.067 0.105 0.057 0.066 0.157 Table 4.6: Percentages of rejections for testing symmetry, with 1000 simulations. 4.4. CRACKS IN CEMENTED FEMORAL COMPONENTS 105 λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 k0g0=fVM1 Q(n);µ;g0 31 0.044 0.072 0.056 0.048 0.047 0.049 0.041 0.058 0.092 0.050 0.081 0.153 Q∗(n);µ 30.043 0.070 0.056 0.043 0.048 0.049 0.046 0.061 0.091 0.053 0.083 0.150 Q(n);g0 30.045 0.069 0.064 0.046 0.040 0.038 0.051 0.046 0.056 0.053 0.056 0.088 Q∗(n) fC0.25 ;3 0.025 0.040 0.062 0.028 0.034 0.033 0.024 0.037 0.059 0.027 0.048 0.094 Q(n);µ;g0 32 0.044 0.072 0.056 0.056 0.096 0.298 0.088 0.240 0.818 0.148 0.480 0.994 Q∗(n);µ 30.043 0.070 0.056 0.065 0.096 0.299 0.090 0.245 0.818 0.149 0.474 0.994 Q(n);g0 30.045 0.069 0.064 0.054 0.074 0.220 0.059 0.146 0.530 0.080 0.221 0.672 Q∗(n) fC0.25 ;3 0.025 0.040 0.062 0.029 0.052 0.247 0.035 0.159 0.697 0.074 0.282 0.924 Q(n);µ;g0 33 0.044 0.072 0.056 0.118 0.298 0.897 0.313 0.812 1 0.651 0.990 1 Q∗(n);µ 30.043 0.070 0.056 0.124 0.299 0.896 0.318 0.810 1 0.660 0.992 1 Q(n);g0 30.045 0.069 0.064 0.091 0.231 0.875 0.196 0.700 1 0.397 0.936 1 Q∗(n) fC0.25 ;3 0.025 0.040 0.062 0.047 0.180 0.850 0.096 0.581 1 0.226 0.797 1 k0g0=fVM10 Q(n);µ;g0 31 0.054 0.046 0.040 0.055 0.084 0.255 0.081 0.209 0.746 0.138 0.440 0.971 Q∗(n);µ 30.051 0.046 0.041 0.055 0.088 0.255 0.082 0.214 0.746 0.135 0.432 0.972 Q(n);g0 30.044 0.051 0.033 0.040 0.048 0.035 0.043 0.043 0.037 0.043 0.038 0.044 Q∗(n) fC0.25 ;3 0.028 0.045 0.040 0.029 0.037 0.033 0.029 0.033 0.038 0.028 0.042 0.053 Q(n);µ;g0 32 0.054 0.046 0.040 0.081 0.177 0.680 0.206 0.553 0.997 0.394 0.902 1 Q∗(n);µ 30.051 0.046 0.041 0.082 0.182 0.682 0.213 0.563 0.997 0.412 0.905 1 Q(n);g0 30.044 0.051 0.033 0.050 0.052 0.043 0.044 0.046 0.066 0.035 0.040 0.049 Q∗(n) fC0.25 ;3 0.028 0.045 0.040 0.030 0.045 0.044 0.041 0.042 0.064 0.041 0.049 0.066 Q(n);µ;g0 33 0.054 0.046 0.040 0.102 0.247 0.858 0.287 0.758 0.999 0.568 0.985 1 Q∗(n);µ 30.051 0.046 0.041 0.101 0.250 0.859 0.289 0.758 0.999 0.590 0.985 1 Q(n);g0 30.044 0.051 0.033 0.057 0.060 0.143 0.049 0.106 0.399 0.061 0.139 0.513 Q∗(n) fC0.25 ;3 0.028 0.045 0.040 0.037 0.057 0.152 0.035 0.085 0.408 0.037 0.095 0.512 k0g0=fC0.45 Q(n);µ;g0 31 0.056 0.052 0.046 0.059 0.044 0.050 0.053 0.040 0.044 0.049 0.043 0.041 Q∗(n);µ 30.052 0.052 0.047 0.054 0.044 0.050 0.049 0.042 0.042 0.049 0.045 0.041 Q(n);g0 30.047 0.053 0.043 0.050 0.041 0.057 0.046 0.035 0.047 0.048 0.032 0.035 Q∗(n) fC0.25 ;3 0.031 0.041 0.039 0.028 0.036 0.054 0.026 0.029 0.040 0.030 0.021 0.033 Q(n);µ;g0 32 0.056 0.052 0.046 0.061 0.092 0.314 0.108 0.260 0.819 0.186 0.473 0.992 Q∗(n);µ 30.052 0.052 0.047 0.061 0.092 0.314 0.118 0.250 0.819 0.174 0.476 0.992 Q(n);g0 30.047 0.053 0.043 0.061 0.075 0.268 0.098 0.181 0.645 0.119 0.247 0.754 Q∗(n) fC0.25 ;3 0.031 0.041 0.039 0.030 0.064 0.283 0.073 0.187 0.786 0.098 0.337 0.962 Q(n);µ;g0 33 0.056 0.052 0.046 0.134 0.307 0.897 0.361 0.814 1 0.660 0.993 1 Q∗(n);µ 30.052 0.052 0.047 0.140 0.306 0.898 0.363 0.811 1 0.674 0.993 1 Q(n);g0 30.047 0.053 0.043 0.100 0.246 0.889 0.224 0.694 1 0.429 0.938 1 Q∗(n) fC0.25 ;3 0.031 0.041 0.039 0.054 0.195 0.889 0.127 0.586 1 0.260 0.838 1 Table 4.7: Percentages of rejections for testing symmetry, with 1000 simulations. 112 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 k0g0=fWN0.5 Q∗(n) fC0.45 ;2 1 0.028 0.045 0.054 0.040 0.060 0.161 0.035 0.081 0.289 0.037 0.067 0.207 Q∗(n) fWC0.5;2 0.027 0.045 0.055 0.039 0.052 0.155 0.037 0.083 0.277 0.038 0.062 0.190 Q∗(n);ˆµ(n) 20.021 0.033 0.032 0.028 0.034 0.107 0.016 0.045 0.152 0.012 0.024 0.069 Q∗(n) fC0.45 ;2 2 0.028 0.045 0.054 0.077 0.253 0.870 0.186 0.715 1 0.363 0.927 1 Q∗(n) fWC0.5;2 0.027 0.045 0.055 0.077 0.247 0.866 0.185 0.724 1 0.344 0.941 1 Q∗(n);ˆµ(n) 20.021 0.033 0.032 0.056 0.174 0.796 0.153 0.643 1 0.324 0.908 1 Q∗(n) fC0.45 ;2 3 0.028 0.045 0.054 0.063 0.103 0.321 0.076 0.213 0.858 0.135 0.451 0.997 Q∗(n) fWC0.5;2 0.027 0.045 0.055 0.054 0.097 0.305 0.072 0.201 0.829 0.125 0.434 0.994 Q∗(n);ˆµ(n) 20.021 0.033 0.032 0.039 0.076 0.317 0.070 0.223 0.886 0.131 0.491 1 k0g0=fWN0.9 Q∗(n) fC0.45 ;2 1 0.024 0.046 0.053 0.031 0.042 0.054 0.026 0.038 0.061 0.032 0.035 0.063 Q∗(n) fWC0.5;2 0.027 0.042 0.054 0.037 0.038 0.063 0.032 0.040 0.043 0.033 0.047 0.050 Q∗(n);ˆµ(n) 200 0 00 0 00 0 00 0 Q∗(n) fC0.45 ;2 2 0.024 0.046 0.053 0.041 0.049 0.127 0.030 0.059 0.360 0.022 0.095 0.467 Q∗(n) fWC0.5;2 0.027 0.042 0.054 0.036 0.041 0.061 0.031 0.043 0.072 0.033 0.038 0.045 Q∗(n);ˆµ(n) 200 0 00 0 00 0 00 0 Q∗(n) fC0.45 ;2 3 0.024 0.046 0.053 0.047 0.113 0.481 0.076 0.287 0.952 0.117 0.485 0.996 Q∗(n) fWC0.5;2 0.027 0.042 0.054 0.041 0.059 0.170 0.047 0.094 0.280 0.059 0.078 0.169 Q∗(n);ˆµ(n) 200 0 00 0 00 0 00 0.028 k0g0=fWC0.5 Q∗(n) fC0.45 ;2 1 0.038 0.058 0.045 0.053 0.070 0.203 0.065 0.152 0.606 0.107 0.347 0.945 Q∗(n) fWC0.5;2 0.040 0.053 0.046 0.046 0.066 0.167 0.060 0.133 0.534 0.091 0.301 0.905 Q∗(n);ˆµ(n) 20.041 0.056 0.046 0.040 0.064 0.218 0.056 0.149 0.676 0.089 0.325 0.967 Q∗(n) fC0.45 ;2 2 0.038 0.058 0.045 0.058 0.129 0.431 0.103 0.303 0.902 0.169 0.469 0.986 Q∗(n) fWC0.5;2 0.040 0.053 0.046 0.056 0.132 0.468 0.102 0.343 0.953 0.169 0.558 0.997 Q∗(n);ˆµ(n) 20.041 0.056 0.046 0.055 0.123 0.432 0.109 0.350 0.932 0.206 0.574 0.995 Q∗(n) fC0.45 ;2 3 0.038 0.058 0.045 0.046 0.062 0.059 0.061 0.069 0.107 0.064 0.064 0.149 Q∗(n) fWC0.5;2 0.040 0.053 0.046 0.044 0.066 0.053 0.055 0.064 0.100 0.047 0.067 0.174 Q∗(n);ˆµ(n) 20.041 0.056 0.046 0.047 0.070 0.154 0.069 0.129 0.449 0.074 0.226 0.762 Table 4.14: Percentages of rejections for testing symmetry, with 1000 simulations. 4.4. CRACKS IN CEMENTED FEMORAL COMPONENTS 113 λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 Test k0g0=fVM1 φ(n);π/32 210.062 0.052 0.101 0.073 0.134 0.540 0.122 0.317 0.933 0.171 0.550 0.999 φ(n);π/16 20.071 0.078 0.228 0.089 0.184 0.732 0.145 0.388 0.968 0.208 0.612 1 φ(n);π/8 20.095 0.178 0.656 0.122 0.298 0.934 0.177 0.486 0.988 0.239 0.646 1 φ(n);π/4 20.126 0.336 0.911 0.126 0.290 0.927 0.120 0.312 0.925 0.126 0.325 0.920 φ(n);π/32 220.062 0.052 0.101 0.142 0.394 0.966 0.388 0.883 1 0.694 0.996 1 φ(n);π/16 20.071 0.078 0.228 0.162 0.471 0.988 0.404 0.892 1 0.706 0.997 1 φ(n);π/8 20.095 0.178 0.656 0.191 0.531 0.997 0.375 0.866 1 0.600 0.990 1 φ(n);π/4 20.126 0.336 0.911 0.126 0.310 0.931 0.129 0.318 0.917 0.143 0.316 0.932 φ(n);π/32 230.062 0.052 0.101 0.074 0.152 0.539 0.132 0.342 0.931 0.197 0.560 0.996 φ(n);π/16 20.071 0.078 0.228 0.090 0.216 0.730 0.151 0.398 0.977 0.208 0.627 0.997 φ(n);π/8 20.095 0.178 0.656 0.124 0.321 0.932 0.168 0.453 0.993 0.226 0.644 0.998 φ(n);π/4 20.126 0.336 0.911 0.122 0.314 0.915 0.131 0.330 0.917 0.142 0.335 0.929 Test k0g0=fVM10 φ(n);π/32 210.355 0.820 1 0.475 0.960 1 0.607 0.991 1 0.732 0.998 1 φ(n);π/16 20.908 1 1 0.942 1 1 0.979 1 1 0.989 1 1 φ(n);π/8 211 1 11 1 11 1 11 1 φ(n);π/4 211 1 11 1 11 1 11 1 φ(n);π/32 220.355 0.820 1 0.580 0.980 1 0.792 0.998 1 0.929 1 1 φ(n);π/16 20.908 1 1 0.967 1 1 0.991 1 1 11 1 φ(n);π/8 211 1 11 1 11 1 11 1 φ(n);π/4 211 1 11 1 11 1 11 1 φ(n);π/32 230.355 0.820 1 0.623 0.989 1 0.856 1 1 0.965 1 1 φ(n);π/16 20.908 1 1 0.976 1 1 0.998 1 1 11 1 φ(n);π/8 211 1 11 1 11 1 11 1 φ(n);π/4 211 1 11 1 11 1 11 1 Test k0g0=fC0.45 φ(n);π/32 210.058 0.057 0.041 0.066 0.089 0.277 0.106 0.254 0.803 0.171 0.475 0.991 φ(n);π/16 20.051 0.052 0.042 0.063 0.084 0.248 0.096 0.228 0.748 0.151 0.430 0.983 φ(n);π/8 20.044 0.055 0.043 0.060 0.071 0.160 0.088 0.152 0.502 0.117 0.269 0.854 φ(n);π/4 20.051 0.054 0.042 0.057 0.048 0.042 0.047 0.053 0.051 0.051 0.045 0.048 φ(n);π/32 220.058 0.057 0.041 0.108 0.283 0.894 0.325 0.810 1 0.641 0.987 1 φ(n);π/16 20.051 0.052 0.042 0.110 0.248 0.861 0.279 0.765 1 0.596 0.980 1 φ(n);π/8 20.044 0.055 0.043 0.079 0.166 0.611 0.186 0.528 0.997 0.378 0.868 1 φ(n);π/4 20.051 0.054 0.042 0.049 0.047 0.046 0.047 0.054 0.060 0.048 0.050 0.046 φ(n);π/32 230.058 0.057 0.041 0.063 0.092 0.285 0.107 0.229 0.807 0.150 0.465 0.987 φ(n);π/16 20.051 0.052 0.042 0.062 0.088 0.255 0.094 0.204 0.736 0.133 0.412 0.975 φ(n);π/8 20.044 0.055 0.043 0.054 0.069 0.149 0.072 0.134 0.464 0.090 0.241 0.815 φ(n);π/4 20.051 0.054 0.042 0.039 0.043 0.046 0.049 0.049 0.055 0.045 0.051 0.052 Table 4.15: Percentages of rejections for testing symmetry, with 1000 simulations. 114 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA λ0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 Test k0g0=fWN0.5 φ(n);π/32 210.058 0.055 0.064 0.071 0.151 0.496 0.125 0.343 0.940 0.231 0.638 0.999 φ(n);π/16 20.058 0.062 0.112 0.070 0.177 0.618 0.121 0.372 0.963 0.219 0.643 1 φ(n);π/8 20.073 0.095 0.267 0.075 0.199 0.704 0.126 0.352 0.959 0.195 0.571 0.999 φ(n);π/4 20.075 0.137 0.493 0.076 0.119 0.510 0.070 0.122 0.503 0.076 0.120 0.482 φ(n);π/32 220.058 0.055 0.064 0.122 0.352 0.939 0.334 0.857 1 0.677 0.996 1 φ(n);π/16 20.058 0.062 0.112 0.133 0.377 0.960 0.327 0.848 1 0.655 0.995 1 φ(n);π/8 20.073 0.095 0.267 0.125 0.367 0.960 0.264 0.739 1 0.488 0.971 1 φ(n);π/4 20.075 0.137 0.493 0.071 0.134 0.488 0.078 0.123 0.485 0.072 0.126 0.489 φ(n);π/32 230.058 0.055 0.064 0.068 0.146 0.483 0.126 0.317 0.923 0.207 0.590 1 φ(n);π/16 20.058 0.062 0.112 0.076 0.168 0.599 0.126 0.334 0.946 0.219 0.567 1 φ(n);π/8 20.073 0.095 0.267 0.077 0.193 0.688 0.118 0.339 0.935 0.185 0.489 0.997 φ(n);π/4 20.075 0.137 0.493 0.073 0.132 0.493 0.070 0.131 0.491 0.066 0.127 0.492 Test k0g0=fWN0.9 φ(n);π/32 210.184 0.490 0.992 0.299 0.820 1 0.450 0.963 1 0.652 0.996 1 φ(n);π/16 20.552 0.982 1 0.724 1 1 0.855 1 1 0.943 1 1 φ(n);π/8 20.990 1 1 0.993 1 1 0.997 1 1 11 1 φ(n);π/4 211 1 11 1 11 1 11 1 φ(n);π/32 220.184 0.490 0.992 0.402 0.920 1 0.680 0.997 1 0.913 1 1 φ(n);π/16 20.552 0.982 1 0.796 0.999 1 0.951 1 1 0.995 1 1 φ(n);π/8 20.990 1 1 0.994 1 1 11 1 11 1 φ(n);π/4 211 1 11 1 11 1 11 1 φ(n);π/32 230.184 0.490 0.992 0.398 0.928 1 0.674 0.996 1 0.903 1 1 φ(n);π/16 20.552 0.982 1 0.782 1 1 0.938 1 1 0.991 1 1 φ(n);π/8 20.990 1 1 0.997 1 1 11 1 11 1 φ(n);π/4 211 1 11 1 11 1 11 1 Test k0g0=fWC0.5 φ(n);π/32 210.076 0.110 0.349 0.103 0.253 0.780 0.153 0.439 0.973 0.221 0.648 1 φ(n);π/16 20.133 0.280 0.894 0.171 0.473 0.986 0.243 0.678 1 0.337 0.838 1 φ(n);π/8 20.295 0.739 1 0.342 0.833 1 0.425 0.912 1 0.508 0.966 1 φ(n);π/4 20.493 0.959 1 0.478 0.944 1 0.471 0.946 1 0.483 0.953 1 φ(n);π/32 220.076 0.110 0.349 0.203 0.539 0.997 0.461 0.943 1 0.771 1 1 φ(n);π/16 20.133 0.280 0.894 0.312 0.745 1 0.569 0.983 1 0.843 1 1 φ(n);π/8 20.295 0.739 1 0.480 0.926 1 0.684 0.996 1 0.863 1 1 φ(n);π/4 20.493 0.959 1 0.480 0.951 1 0.486 0.951 1 0.484 0.954 1 φ(n);π/32 230.076 0.110 0.349 0.107 0.308 0.851 0.190 0.537 0.994 0.300 0.756 1 φ(n);π/16 20.133 0.280 0.894 0.181 0.518 0.996 0.276 0.738 1 0.396 0.880 1 φ(n);π/8 20.295 0.739 1 0.353 0.828 1 0.422 0.923 1 0.520 0.964 1 φ(n);π/4 20.493 0.959 1 0.485 0.942 1 0.472 0.948 1 0.473 0.946 1 Table 4.16: Percentages of rejections for testing symmetry, with 1000 simulations. 4.4. CRACKS IN CEMENTED FEMORAL COMPONENTS 115 ν0 0.2 0.4 0.6 n30 100 500 30 100 500 30 100 500 30 100 500 KJ10(0,0.5,0.5,ν) Q∗(n);µ 20.053 0.051 0.058 0.194 0.538 0.995 0.524 0.962 1 0.747 0.997 1 Q∗(n) fVMκ;2 0.044 0.057 0.050 0.047 0.057 0.040 0.053 0.059 0.047 0.044 0.061 0.067 KJ10(0,0.9,0.5,ν) Q∗(n);µ 20.048 0.038 0.059 0.271 0.724 1 0.707 0.997 1 0.916 1 1 Q∗(n) fVMκ;2 0.049 0.042 0.053 0.053 0.048 0.062 0.052 0.067 0.085 0.064 0.070 0.142 KJ15(0,0.5,ν) Q∗(n);µ 20.054 0.043 0.047 0.070 0.103 0.399 0.127 0.307 0.923 0.199 0.569 0.999 Q∗(n) fVMκ;2 0.038 0.053 0.053 0.069 0.113 0.364 0.103 0.296 0.895 0.175 0.567 1 KJ15(0,0.9,ν) Q∗(n);µ 20.064 0.049 0.036 0.051 0.067 0.165 0.066 0.154 0.496 0.123 0.263 0.791 Q∗(n) fVMκ;2 0.025 0.045 0.062 0.032 0.058 0.173 0.048 0.157 0.516 0.061 0.279 0.867 Table 4.17: Percentages of rejections for testing symmetry, with 1000 simulations. for that reason, it was also excluded from their study, which aimed to check if the cracks were uniformly distributed. This leaves a total number of cement cracks equal to 2388: 1567 in the proximal region, and 434 in the distal region. Circular data plots for the two regions, together with rose diagrams and kernel density estimates obtained using the plug–in rule of Oliveira et al. (2012) to select the concentration parameter, are portrayed in Figure 4.2. Since neither the mean direction, nor the underlying distribution, nor the value of kare assumed to be known, the φ∗(n) fVMκ;2 test was applied as it is a powerful omnibus test under such circumstances. The p–value obtained for the directions of the cracks in the proximal region was 0.6096, whilst for those in the distal region the p–value was 0.0483. Clearly, there is no significant evidence to reject underlying reflective symmetry for the distribution of the directions of the cracks in the proximal region. On the other hand, for a significance level of the 5% the test rejects reflective symmetry for the distribution of the crack directions in the distal region. Observing the p–values obtained, for other values of k, when the φ∗(n) fC0.25 ;ktest is used, in the proximal region there are still no evidences for rejecting the reflective symmetry, as for k= 1 the p–value is 0.7543, for k= 2 it is 0.6058 and for k= 3 it is 0.5730. Now, in the distal region, the conclusion is not so clear as the p–values are 0.0345 (for k= 1), 0.0451 (for k= 2) and 0.9507 (for k= 3). The distal sample can be seen as another good example of why the multimodality test are not only important per se, applying the new test proposed in Section 3.2 for H0:j= 1 (with B= 1000), the p–value was 0.259 so there is no evidence for rejecting the null hypothesis of unimodality. As, usually, the k–sine–skewed are unimodal when k= 1 (if the data is not “too” concentrated 116 CHAPTER 4. SYMMETRY TESTING IN CIRCULAR DATA around the central direction), in this case, it seems that the most powerful test can be obtained with a low value of k, which supports the evidences against the null for a significance level of the 5%. These findings provide further evidence to complement the results presented in Mann et al. (2003) regarding the different distributions of the crack directions in the two regions. It seems that cracks in the proximal and distal region are not equally distributed just with a translation in the central direction. Chapter 5 Conclusions and discussion This essay was focused on assessing simplifying hypothesis in density estimation for linear and circular data. The main statistical contributions that have been achieved in this PhD thesis are summed up in the following points. Testing the number of modes for linear data. The objective was to provide a nonparametric tool for testing if the linear random variable has kmodes. In addition, the proposed procedure has been compared with the existing methods, outperforming all of them in a wide range of scenarios. This goal corresponds to the contents in Chapter 2. Testing the number of modes for circular data. The second goal was to extend the previous testing procedure to the circular setting. Since the final objective of this part was to determine if, in each region of the world, there is one season of fires or more, it was also studied how to handle the multiple testing problem with spatially correlated data. The results achieved in this point are collected in Chapter 3. Testing reflective symmetry for circular data. The third objective was to introduce a new proposal for testing circular symmetry when the central location is unknown. The development of this point is provided in Chapter 4. In what follows, some final comments and discussion on each chapter are presented. 117 118 CHAPTER 5. CONCLUSIONS AND DISCUSSION 5.1 Multimodality tests for linear data Determining the number of modes in a distribution is a relevant practical problem in many applied sciences. The proposal presented in Chapter 2 provides a good performance for the testing problem (1.1.3), being in the case of a general number of modes kthe only alternative with a reasonable behaviour. The totally nonparametric testing procedure can be extended to other contexts where a natural nonparametric estimator under the null hypothesis is available, as it was done in Chapter 3 for circular data. The relevance of the developed techniques was shown in two datasets. First, in a classical example in the philatelic field, given an answer about the number of modes with a well–calibrated test. Second, showing the method potential for determining if the malaria is not eradicated in a given population. With the aim of making this procedure accessible for the scientific community, and therefore, enabling its use in other practical problems, the Rpackage multimode was developed and it is presented in Appendix C. Future research on this topic includes the extension of this techniques to the d– dimensional case (with d > 1). The excess mass statistic can be adapted in this setting, but there is not a clear definition of the critical bandwidth. In fact, using the Gaussian kernel, for the simplest case in the bidimensional setting, where a single value of his taken in both dimensions, in Scott (2015, Sect. 9.2.4), an counterexample where the monotony in the number of modes is not satisfied can be found. Then, empirical approximations must be developed for tackling this problem. 5.2 Multimodality tests for circular data As mentioned in the Introduction, asymmetric and multimodal distribution occur frequently in practice in the circular setting. In Chapters 3 and 4 different tests were developed for testing both hypothesis in a flexible way. The growing interest in the last few years in more flexible models in circular data, give to these methods a special relevance as a preliminary tool before applying more complicated methods. The proposal presented in Chapter 3 shows how the method presented in Chapter 2 can be extended to other settings, having a correct behaviour in the circular setting. The detailed method, adapting the proposal of Benjamini and Heller (2007), provides a useful algorithm for correcting the FDR accounting for the spatial dependence of the data in other contexts where prior information about the neighbouring 5.3. SYMMETRY TESTS IN THE CIRCULAR SETTING 119 locations is known. Related with the last point, future research on this topic includes trying to integrate the spatial correlation in the test statistic before making the FDR correction. The map shown in Section 3.4.3 provides new conclusions about how human activity modifies fire seasonality. Related with this application, as it is expected that the circular kernel density estimator provides good estimation of the location of the modes and antimodes, future research may include the use of multimodality test as a preliminary tool for exploring when the peaks of fires are produced and also its associated mass. This will allow review different works in the forest field with nonparametric techniques. For instance, one can determine when the principal peaks of fires are produced in each 0.5◦cell (Le Page et al., 2010), the delay of the agricultural fires with respect to the climatological ones (Magi et al., 2012) or the mass associated to each peak for better understanding the importance of the different human activities (Korontzi et al., 2006). Other possible approach for solving the last questions can be done using a mixture of flexible parametric distributions. In particular, a future objective includes trying to model the wildfires distribution using Hidden Markov Models (HMM). The conclusions about the number of modes can be used as a preliminary tool for determining how many components should be used in the mixture distribution used in the HMM. Finally, in future updates of the Rpackage, the objective will be include different tools for determining the number of modes in the circular setting. 5.3 Symmetry tests in the circular setting Test for circular reflective symmetry have been developed in Chapter 4. Specifically, such tests consider an unknown centre of symmetry and are optimal against k–sine– skewed alternatives. Recommendations for their use, as well as other tests that have been proposed in the literature, were established in the light of the simulation based results reported in Section 4.3. As mentioned there, the proposed tests are generally conservative when the sample size is of the order of 30. For that reason, future research may involve more sophisticated methods for generating the bootstrap resamples in order to better approximate the distribution of the test statistic under the null (improving the bootstrap method mentioned in Section 4.3) and also for having a complete non–parametric approach. Circular data are just one class of directional data. Others include bivariate circular data distributed on the torus, cylindrical data, spherical data and data distributed 120 CHAPTER 5. CONCLUSIONS AND DISCUSSION on the surfaces of the extensions of such Riemannian manifolds. The development of tests for reflective symmetry on such manifolds would be of considerable interest. Ideas underpinning such tests are explored in Jupp and Spurr (1983) and Jupp et al. (2016). Appendices 121 128 APPENDIX A. MODELS FOR SIMULATION STUDIES −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M1 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M5 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 x f(x) M9 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M2 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M6 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M10 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) M3 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) M7 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) M26 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M4 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M8 Figure A.1: Unimodal linear density functions: M1–M10 and M26. These models were used in Section 1.1 for illustration purposes and in Section 2.2 as part of the simulation study carried out for comparing the calibration of the different procedures for testing the number of modes in the linear case. 129 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M11 −0.5 0.0 0.5 1.0 1.5 0 1 2 3 x f(x) M15 −0.5 0.0 0.5 1.0 1.5 0 1 2 3 4 5 6 x f(x) M19 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M23 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 x f(x) M12 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 2.5 x f(x) M16 −0.5 0.0 0.5 1.0 1.5 0 1 2 3 4 5 6 x f(x) M20 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M24 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 2.5 x f(x) M13 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M17 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M21 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) M25 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) M14 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 x f(x) M18 −0.5 0.0 0.5 1.0 1.5 0.0 0.5 1.0 1.5 2.0 x f(x) M22 Figure A.2: Linear density functions. M11–M20: bimodal models. M21–M25: trimodal models. These models were used in Sections 1.1 and 2.1.3 for illustration purposes and in Section 2.2 as part of the simulation study comparing, in terms of empirical size and power, the different procedures for testing the number of modes in the linear case. 130 APPENDIX A. MODELS FOR SIMULATION STUDIES 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 θ f(θ) 0π2π3π22π MC1 0.0 0.2 0.4 0.6 0.8 1.0 1.2 θ f(θ) 0π2π3π22π MC5 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 θ f(θ) MC9 0.0 0.2 0.4 0.6 0.8 θ f(θ) 0π2π3π22π MC2 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π MC6 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 θ f(θ) MC10 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 θ f(θ) 0π2π3π22π MC3 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π MC7 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π MC4 0.00 0.05 0.10 0.15 0.20 0.25 0.30 θ f(θ) 0π2π3π22π MC8 Figure A.3: Linear representation of the unimodal circular denstity functions: MC1– MC10. These models were used in Section 1.2 and 3.2 for illustration purposes and in Section 3.3 as part of the simulation study carried out for comparing the calibration of the different procedures for testing the number of modes in the circular case. 131 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π MC11 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) MC15 0.00 0.05 0.10 0.15 0.20 θ f(θ) 0π2π3π22π MC19 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 θ f(θ) 0π2π3π22π MC23 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π MC12 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) MC16 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 θ f(θ) 0π2π3π22π MC20 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π MC24 0.00 0.05 0.10 0.15 0.20 0.25 θ f(θ) 0π2π3π22π MC13 0.0 0.1 0.2 0.3 θ f(θ) 0π2π3π22π MC17 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π MC21 0 1 2 3 4 5 6 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) MC25 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π MC14 0.0 0.1 0.2 0.3 0.4 θ f(θ) 0π2π3π22π MC18 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) 0π2π3π22π MC22 Figure A.4: Linear representation of circular density functions. MC11–MC20: bimodal models. MC21–MC25: trimodal models. These models were used in Section 3.3 as part of the simulation study comparing, in terms of empirical size and power, the different procedures for testing the number of modes in the circular case. 132 APPENDIX A. MODELS FOR SIMULATION STUDIES 0 π 2 π 3π 2 + MC1 0 π 2 π 3π 2 + MC5 0 π 2 π 3π 2 + MC9 0 π 2 π 3π 2 + MC2 0 π 2 π 3π 2 + MC6 0 π 2 π 3π 2 + MC10 0 π 2 π 3π 2 + MC3 0 π 2 π 3π 2 + MC7 0 π 2 π 3π 2 + MC4 0 π 2 π 3π 2 + MC8 Figure A.5: Circular representation of the unimodal density functions: MC1–MC10. These models were used in Section 1.2 and 3.2 for illustration purposes and in Section 3.3 as part of the simulation study carried out for comparing the calibration of the different procedures for testing the number of modes in the circular case. 133 0 π 2 π 3π 2 + MC11 0 π 2 π 3π 2 + MC15 0 π 2 π 3π 2 + MC19 0 π 2 π 3π 2 + MC23 0 π 2 π 3π 2 + MC12 0 π 2 π 3π 2 + MC16 0 π 2 π 3π 2 + MC20 0 π 2 π 3π 2 + MC24 0 π 2 π 3π 2 + MC13 0 π 2 π 3π 2 + MC17 0 π 2 π 3π 2 + MC21 0 π 2 π 3π 2 + MC25 0 π 2 π 3π 2 + MC14 0 π 2 π 3π 2 + MC18 0 π 2 π 3π 2 + MC22 Figure A.6: Circular representation of some density functions. MC11–MC20: bimodal models. MC21–MC25: trimodal models. These models were used in Section 3.3 as part of the simulation study comparing, in terms of empirical size and power, the different procedures for testing the number of modes in the circular case. 134 APPENDIX A. MODELS FOR SIMULATION STUDIES 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssVM(0,1,0,k0) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 θ f(θ) π3π20π2π kssVM(0,10,0,k0) 0.0 0.1 0.2 0.3 0.4 θ f(θ) π3π20π2π kssC(0,0.45,0,k0) 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssWN(0,0.5,0,k0) 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssVM(0,1,0.2,k0) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 θ f(θ) π3π20π2π kssVM(0,10,0.2,k0) 0.0 0.1 0.2 0.3 0.4 θ f(θ) π3π20π2π kssC(0,0.45,0.2,k0) 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssWN(0,0.5,0.2,k0) 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssVM(0,1,0.4,k0) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 θ f(θ) π3π20π2π kssVM(0,10,0.4,k0) 0.0 0.1 0.2 0.3 0.4 θ f(θ) π3π20π2π kssC(0,0.45,0.4,k0) 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssWN(0,0.5,0.4,k0) 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssVM(0,1,0.6,k0) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 θ f(θ) π3π20π2π kssVM(0,10,0.6,k0) 0.0 0.1 0.2 0.3 0.4 θ f(θ) π3π20π2π kssC(0,0.45,0.6,k0) 0.0 0.1 0.2 0.3 0.4 0.5 θ f(θ) π3π20π2π kssWN(0,0.5,0.6,k0) Figure A.7: Linear representation of k–sine–skewed models. solid line: k0= 1. Dashed line: k0= 2. Pointed line: k0= 3. These models were used in Section 4.3 as part of the simulation study comparing, in terms of empirical size and power, the different procedures for testing reflective symmetry in the circular case. 135 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π kssWN(0,0.9,0,k0) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) π3π20π2π kssWC(0,0.5,0,k0) 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π KJ10(0, κ, 0.5,0) 0.0 0.5 1.0 1.5 2.0 2.5 3.0 θ f(θ) π3π20π2π KJ15(0, γ, 0) 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π kssWN(0,0.9,0.2,k0) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) π3π20π2π kssWC(0,0.5,0.2,k0) 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π KJ10(0, κ, 0.5,0.2) 0.0 0.5 1.0 1.5 2.0 2.5 3.0 θ f(θ) π3π20π2π KJ15(0, γ, 0.2) 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π kssWN(0,0.9,0.4,k0) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) π3π20π2π kssWC(0,0.5,0.4,k0) 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π KJ10(0, κ, 0.5,0.4) 0.0 0.5 1.0 1.5 2.0 2.5 3.0 θ f(θ) π3π20π2π KJ15(0, γ, 0.4) 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π kssWN(0,0.9,0.6,k0) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 θ f(θ) π3π20π2π kssWC(0,0.5,0.6,k0) 0.0 0.2 0.4 0.6 0.8 1.0 θ f(θ) π3π20π2π KJ10(0, κ, 0.5,0.6) 0.0 0.5 1.0 1.5 2.0 2.5 3.0 θ f(θ) π3π20π2π KJ15(0, γ, 0.6) Figure A.8: Linear representation of k–sine–skewed and Kato Jones models. First two columns: solid line, k0= 1; dashed line, k0= 2; pointed line, k0= 3. Third column: solid line, κ= 0.5; dashed line, κ= 0.9. Forth column: solid line, γ= 0.5; dashed line, γ= 0.9. These models were used in Section 4.3 as part of the simulation study comparing, in terms of empirical size and power, the different procedures for testing reflective symmetry in the circular case. 136 APPENDIX A. MODELS FOR SIMULATION STUDIES 0 π 2 π 3π 2 + kssVM(0,1,0,k0) 0 π 2 π 3π 2 + kssVM(0,10,0,k0) 0 π 2 π 3π 2 + kssC(0,0.45,0,k0) 0 π 2 π 3π 2 + kssWN(0,0.5,0,k0) 0 π 2 π 3π 2 + kssVM(0,1,0.2,k0) 0 π 2 π 3π 2 + kssVM(0,10,0.2,k0) 0 π 2 π 3π 2 + kssC(0,0.45,0.2,k0) 0 π 2 π 3π 2 + kssWN(0,0.5,0.2,k0) 0 π 2 π 3π 2 + kssVM(0,1,0.4,k0) 0 π 2 π 3π 2 + kssVM(0,10,0.4,k0) 0 π 2 π 3π 2 + kssC(0,0.45,0.4,k0) 0 π 2 π 3π 2 + kssWN(0,0.5,0.4,k0) 0 π 2 π 3π 2 + kssVM(0,1,0.6,k0) 0 π 2 π 3π 2 + kssVM(0,10,0.6,k0) 0 π 2 π 3π 2 + kssC(0,0.45,0.6,k0) 0 π 2 π 3π 2 + kssWN(0,0.5,0.6,k0) Figure A.9: Circular representation of k–sine–skewed models. solid line: k0= 1. Dashed line: k0= 2. Pointed line: k0= 3. These models were used in Section 4.3 as part of the simulation study comparing, in terms of empirical size and power, the different procedures for testing reflective symmetry in the circular case. 137 0 π 2 π 3π 2 + kssWN(0,0.9,0,k0) 0 π 2 π 3π 2 + kssWC(0,0.5,0,k0) 0 π 2 π 3π 2 + KJ10(0, κ, 0.5,0) 0 π 2 π 3π 2 + KJ15(0, γ, 0) 0 π 2 π 3π 2 + kssWN(0,0.9,0.2,k0) 0 π 2 π 3π 2 + kssWC(0,0.5,0.2,k0) 0 π 2 π 3π 2 + KJ10(0, κ, 0.5,0.2) 0 π 2 π 3π 2 + KJ15(0, γ, 0.2) 0 π 2 π 3π 2 + kssWN(0,0.9,0.4,k0) 0 π 2 π 3π 2 + kssWC(0,0.5,0.4,k0) 0 π 2 π 3π 2 + KJ10(0, κ, 0.5,0.4) 0 π 2 π 3π 2 + KJ15(0, γ, 0.4) 0 π 2 π 3π 2 + kssWN(0,0.9,0.6,k0) 0 π 2 π 3π 2 + kssWC(0,0.5,0.6,k0) 0 π 2 π 3π 2 + KJ10(0, κ, 0.5,0.6) 0 π 2 π 3π 2 + KJ15(0, γ, 0.6) Figure A.10: Circular representation of k–sine–skewed and Kato Jones models. First two columns: solid line, k0= 1; dashed line, k0= 2; pointed line, k0= 3. Third column: solid line, κ= 0.5; dashed line, κ= 0.9. Forth column: solid line, γ= 0.5; dashed line, γ= 0.9. These models were used in Section 4.3 as part of the simulation study comparing, in terms of empirical size and power, the different procedures for testing reflective symmetry in the circular case.