1 /* 2 * This file is part of ELKI: 3 * Environment for Developing KDD-Applications Supported by Index-Structures 4 * 5 * Copyright (C) 2018 6 * ELKI Development Team 7 * 8 * This program is free software: you can redistribute it and/or modify 9 * it under the terms of the GNU Affero General Public License as published by 10 * the Free Software Foundation, either version 3 of the License, or 11 * (at your option) any later version. 12 * 13 * This program is distributed in the hope that it will be useful, 14 * but WITHOUT ANY WARRANTY; without even the implied warranty of 15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 16 * GNU Affero General Public License for more details. 17 * 18 * You should have received a copy of the GNU Affero General Public License 19 * along with this program. If not, see <http://www.gnu.org/licenses/>. 20 */ 21 package de.lmu.ifi.dbs.elki.math.statistics.distribution.estimator; 22 23 import de.lmu.ifi.dbs.elki.math.MathUtil; 24 import de.lmu.ifi.dbs.elki.math.statistics.distribution.NormalDistribution; 25 import de.lmu.ifi.dbs.elki.math.statistics.distribution.SkewGeneralizedNormalDistribution; 26 import de.lmu.ifi.dbs.elki.utilities.documentation.Reference; 27 import de.lmu.ifi.dbs.elki.utilities.optionhandling.AbstractParameterizer; 28 import net.jafama.FastMath; 29 30 /** 31 * Estimate the parameters of a skew Normal Distribution (Hoskin's Generalized 32 * Normal Distribution), using the methods of L-Moments (LMM). 33 * <p> 34 * Reference: 35 * <p> 36 * J. R. M. Hosking<br> 37 * Fortran routines for use with the method of L-moments Version 3.03<br> 38 * IBM Research. 39 * 40 * @author Erich Schubert 41 * @since 0.6.0 42 * 43 * @has - - - SkewGeneralizedNormalDistribution 44 */ 45 @Reference(authors = "J. R. M. Hosking", // 46 title = "Fortran routines for use with the method of L-moments Version 3.03", // 47 booktitle = "IBM Research Technical Report", // 48 bibkey = "tr/ibm/Hosking00") 49 public class SkewGNormalLMMEstimator implements LMMDistributionEstimator<SkewGeneralizedNormalDistribution> { 50 /** 51 * Static instance. 52 */ 53 public static final SkewGNormalLMMEstimator STATIC = new SkewGNormalLMMEstimator(); 54 55 /** Polynomial approximation */ 56 private static final double // 57 A0 = 0.20466534e+01, // 58 A1 = -0.36544371e+01, // 59 A2 = 0.18396733e+01, // 60 A3 = -0.20360244; 61 62 /** Polynomial approximation */ 63 private static final double // 64 B1 = -0.20182173e+01, // 65 B2 = 0.12420401e+01, // 66 B3 = -0.21741801; 67 68 /** 69 * Constructor. Private: use static instance. 70 */ SkewGNormalLMMEstimator()71 private SkewGNormalLMMEstimator() { 72 super(); 73 } 74 75 @Override getNumMoments()76 public int getNumMoments() { 77 return 3; 78 } 79 80 @Override estimateFromLMoments(double[] xmom)81 public SkewGeneralizedNormalDistribution estimateFromLMoments(double[] xmom) { 82 if(!(xmom[1] > 0.) || !(Math.abs(xmom[2]) < 1.0)) { 83 throw new ArithmeticException("L-Moments invalid"); 84 } 85 // Generalized Normal Distribution estimation: 86 double t3 = xmom[2]; 87 if(Math.abs(t3) >= .95) { 88 // Extreme skewness 89 return new SkewGeneralizedNormalDistribution(0., -1., 0.); 90 } 91 if(Math.abs(t3) <= 1e-8) { 92 // t3 effectively zero, this is a normal distribution. 93 return new SkewGeneralizedNormalDistribution(xmom[0], xmom[1] * MathUtil.SQRTPI, 0.); 94 } 95 final double tt = t3 * t3; 96 double shape = -t3 * (A0 + tt * (A1 + tt * (A2 + tt * A3))) / (1. + tt * (B1 + tt * (B2 + tt * B3))); 97 final double e = FastMath.exp(.5 * shape * shape); 98 double scale = xmom[1] * shape / (e * NormalDistribution.erf(.5 * shape)); 99 double location = xmom[0] + scale * (e - 1.) / shape; 100 return new SkewGeneralizedNormalDistribution(location, scale, shape); 101 } 102 103 @Override getDistributionClass()104 public Class<? super SkewGeneralizedNormalDistribution> getDistributionClass() { 105 return SkewGeneralizedNormalDistribution.class; 106 } 107 108 @Override toString()109 public String toString() { 110 return this.getClass().getSimpleName(); 111 } 112 113 /** 114 * Parameterization class. 115 * 116 * @author Erich Schubert 117 */ 118 public static class Parameterizer extends AbstractParameterizer { 119 @Override makeInstance()120 protected SkewGNormalLMMEstimator makeInstance() { 121 return STATIC; 122 } 123 } 124 } 125