sascha@2645: package de.intevation.flys.artifacts.math;
sascha@2645: 
sascha@2645: import org.apache.commons.math.stat.descriptive.moment.Mean;
sascha@2645: import org.apache.commons.math.stat.descriptive.moment.StandardDeviation;
sascha@2645: 
sascha@2645: import org.apache.commons.math.distribution.TDistributionImpl;
sascha@2645: 
sascha@2645: import java.util.List;
sascha@2645: import java.util.ArrayList;
sascha@2645: 
sascha@2645: public class Outlier
sascha@2645: {
sascha@2645:     public static class IndexedValue {
sascha@2645:         protected int    index;
sascha@2645:         protected double value;
sascha@2645: 
sascha@2645:         public IndexedValue() {
sascha@2645:         }
sascha@2645: 
sascha@2645:         public IndexedValue(int index, double value) {
sascha@2645:             this.index = index;
sascha@2645:             this.value = value;
sascha@2645:         }
sascha@2645: 
sascha@2645:         public int getIndex() {
sascha@2645:             return index;
sascha@2645:         }
sascha@2645: 
sascha@2645:         public void setIndex(int index) {
sascha@2645:             this.index = index;
sascha@2645:         }
sascha@2645: 
sascha@2645:         public double getValue() {
sascha@2645:             return value;
sascha@2645:         }
sascha@2645: 
sascha@2645:         public void setValue(double value) {
sascha@2645:             this.value = value;
sascha@2645:         }
sascha@2645:     } // class IndexedValue
sascha@2645: 
sascha@2645:     public Outlier() {
sascha@2645:     }
sascha@2645: 
sascha@2645:     public static List<IndexedValue> findOutliers(
sascha@2645:         List<IndexedValue> inputValues,
sascha@2645:         double             alpha
sascha@2645:     ) {
sascha@2645:         ArrayList<IndexedValue> outliers = new ArrayList<IndexedValue>();
sascha@2645: 
sascha@2645:         ArrayList<IndexedValue> values =
sascha@2645:             new ArrayList<IndexedValue>(inputValues);
sascha@2645: 
sascha@2645:         for (;;) {
sascha@2645:             int N = values.size();
sascha@2645: 
sascha@2645:             if (N < 4) {
sascha@2645:                 break;
sascha@2645:             }
sascha@2645: 
sascha@2645:             Mean mean = new Mean();
sascha@2645:             StandardDeviation std = new StandardDeviation();
sascha@2645: 
sascha@2645:             for (IndexedValue value: values) {
sascha@2645:                 mean.increment(value.getValue());
sascha@2645:                 std.increment(value.getValue());
sascha@2645:             }
sascha@2645: 
sascha@2645:             double m = mean.getResult();
sascha@2645:             double s = std.getResult();
sascha@2645: 
sascha@2645:             double maxZ = -Double.MAX_VALUE;
sascha@2645:             int iv = -1;
sascha@2645:             for (int i >= 0; i = N-1; --i) {
sascha@2645:                 IndexedValue v = values.get(i);
sascha@2645:                 double z = Math.abs(m - v.getValue())/s;
sascha@2645:                 if (z > maxZ) {
sascha@2645:                     maxZ = z;
sascha@2645:                     iv = i;
sascha@2645:                 }
sascha@2645:             }
sascha@2645: 
sascha@2645:             double t = Math.sqrt((N*(N-2)*z*z)/((N-1)*(N-1) - N*z*z))
sascha@2645: 
sascha@2645:             TDistributionImpl tdist = new TDistributionImpl(N-2);
sascha@2645: 
sascha@2645:             double p = tdist.cumulativeProbability(t)
sascha@2645:         }
sascha@2645: 
sascha@2645: 
sascha@2645:         return outliers;
sascha@2645:     }
sascha@2645: }
sascha@2645: // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :