cern.jet.stat.Descriptive Java Examples

/**
 * Returns the auto-correlation of a data sequence.
 *
 * @param scope
 * @param data
 * @param lag
 * @param mean
 * @param variance
 * @return
 */
@operator (
		value = "auto_correlation",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the auto-correlation of a data sequence given some lag",
		comment = "",
		examples = { @example (
				value = "auto_correlation([1,0,1,0,1,0],2)",
				equals = "1"),
				@example (
						value = "auto_correlation([1,0,1,0,1,0],1)",
						equals = "-1") })
public static Double opAutoCorrelation(final IScope scope, final IContainer data, final Integer lag) {

	// TODO input parameters validation

	final double mean = (Double) Containers.mean(scope, data);
	final double variance = Stats.opVariance(scope, data);

	return Descriptive.autoCorrelation(from(scope, data), lag, mean, variance);
}
 

/**
 * assertion: isBasicParametersValid == false
 * 
 */
protected void updateIncrementalStats() {
	// prepare arguments
	double[] arguments = new double[4];
	arguments[0] = this.min;
	arguments[1] = this.max;
	arguments[2] = this.sum;
	arguments[3] = this.sum_xx;

	Descriptive.incrementalUpdate(this.elements,this.size, this.elements.size()-1,arguments);

	// store the new parameters back
	this.min =   arguments[0];
	this.max =   arguments[1];
	this.sum =  arguments[2];
	this.sum_xx = arguments[3];

	this.isIncrementalStatValid = true;
	this.size = this.elements.size(); // next time we don't need to redo the stuff we have just done...
}
 

/**
 * Adds the part of the specified list between indexes <tt>from</tt> (inclusive) and <tt>to</tt> (inclusive) to the receiver.
 *
 * @param list the list of which elements shall be added.
 * @param from the index of the first element to be added (inclusive).
 * @param to the index of the last element to be added (inclusive).
 * @throws IndexOutOfBoundsException if <tt>list.size()&gt;0 && (from&lt;0 || from&gt;to || to&gt;=list.size())</tt>.
 */
public synchronized void addAllOfFromTo(DoubleArrayList list, int from, int to) {
	super.addAllOfFromTo(list, from, to);
	
	if (this.sumOfPowers != null) {
		//int max_k = this.min_k + this.sumOfPowers.length-1;
		Descriptive.incrementalUpdateSumsOfPowers(list, from, to, 3, getMaxOrderForSumOfPowers(), this.sumOfPowers);
	}

	if (this.hasSumOfInversions) {
		this.sumOfInversions += Descriptive.sumOfInversions(list, from, to);
	}
	
	if (this.hasSumOfLogarithms) {
		this.sumOfLogarithms += Descriptive.sumOfLogarithms(list, from, to);
	}
}
 

/**
 * Adds the part of the specified list between indexes <tt>from</tt> (inclusive) and <tt>to</tt> (inclusive) to the receiver.
 *
 * @param list the list of which elements shall be added.
 * @param from the index of the first element to be added (inclusive).
 * @param to the index of the last element to be added (inclusive).
 * @throws IndexOutOfBoundsException if <tt>list.size()&gt;0 && (from&lt;0 || from&gt;to || to&gt;=list.size())</tt>.
 */
public synchronized void addAllOfFromTo(DoubleArrayList list, int from, int to) {
	//if (this.arguments == null) setUpCache();
	synchronized (arguments) {
		// prepare arguments
		arguments[0] = this.min;
		arguments[1] = this.max;
		arguments[2] = this.sum;
		arguments[3] = this.sum_xx;

		Descriptive.incrementalUpdate(list, from, to, arguments);

		// store the new parameters back
		this.min = arguments[0];
		this.max = arguments[1];
		this.sum = arguments[2];
		this.sum_xx = arguments[3];

		this.size += to-from+1;
	}
}
 

/**
 *
 *
 * @param scope
 * @param size
 * @param sum
 * @param numOfSquares
 * @return
 */
@operator (
		value = "variance",
		can_be_const = true,
		type = IType.FLOAT,
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the variance of a data sequence. That is (sumOfSquares - mean*sum) / size with mean = sum/size.",
		comment = "In the example we consider variance of [1,3,5,7]. The size is 4, the sum is 1+3+5+7=16 and the sum of squares is 84."
				+ "The variance is (84- 16^2/4)/4. CQFD.",
		examples = { @example (
				value = "int(variance(4,16,84))",
				equals = "5",
				returnType = "int") })
public static Double variance(final IScope scope, final Integer size, final Double sum,
		final Double numOfSquares) {

	// TODO input parameters validation

	return Descriptive.variance(size, sum, numOfSquares);
}
 

@operator (
		value = "variance",
		can_be_const = true,
		type = IType.FLOAT,
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the variance from a standard deviation.",
		comment = "",
		examples = { @example (
				value = "int(variance([1,3,5,6,9,11,12,13]))",
				equals = "17",
				returnType = "int") })

// @test ("int(variance([1,3,5,6,9,11,12,13])) = 17")
public static Double opVariance(final IScope scope, final Double standardDeviation) {

	// TODO input parameters validation

	return Descriptive.variance(standardDeviation);
}
 

/**
 *
 *
 * @param scope
 * @param moment3
 * @param standardDeviation
 * @return
 */
@operator (
		value = "skew",
		can_be_const = true,
		type = IType.FLOAT,
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the skew of a data sequence when the 3rd moment has already been computed.",
		comment = "In R moment(c(1, 3, 5, 6, 9, 11, 12, 13), order=3,center=TRUE) is -10.125 and sd(c(1,3,5,6,9,11,12,13)) = 4.407785"
				+ "The value of the skewness tested here is different because there are different types of estimator"
				+ "Joanes and Gill (1998) discuss three methods for estimating skewness:"
				+ "Type 1: g_1 = m_3 / m_2^(3/2). This is the typical definition used in many older textbooks."
				+ "Type 2: G_1 = g_1 * sqrt(n(n-1)) / (n-2). Used in SAS and SPSS."
				+ "Type 3: b_1 = m_3 / s^3 = g_1 ((n-1)/n)^(3/2). Used in MINITAB and BMDP."
				+ "In R skewness(c(1, 3, 5, 6, 9, 11, 12, 13),type=3) is -0.1182316",
		examples = { @example (
				value = "skew(-10.125,4.407785) with_precision(2)",
				equals = "-0.12") })
public static Double opSkew(final IScope scope, final Double moment3, final Double standardDeviation) {

	// TODO input parameters validation

	return Descriptive.skew(moment3, standardDeviation);
}
 

/**
 *
 *
 * @param scope
 * @param data
 * @param mean
 * @param standardDeviation
 * @return
 */
@operator (
		value = "skew",
		can_be_const = true,
		type = IType.FLOAT,
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the skew of a data sequence, which is moment(data,3,mean) / standardDeviation3",
		comment = "",
		examples = { @example (
				value = "skew([1,3,5,6,9,11,12,13]) with_precision(2)",
				equals = "-0.14") })
public static Double opSkew(final IScope scope, final IContainer data) {

	// TODO input parameters validation

	final double mean = (Double) Containers.mean(scope, data);
	final double standardDeviation = Stats.opStDev(scope, data);

	return Descriptive.skew(from(scope, data), mean, standardDeviation);
}
 

/**
 *
 *
 * @param scope
 * @param size
 * @param sumOfSquares
 * @return
 */
@operator (
		value = "rms",
		can_be_const = true,
		type = IType.FLOAT,
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the RMS (Root-Mean-Square) of a data sequence. "
				+ "The RMS of data sequence is the square-root of the mean of the squares "
				+ "of the elements in the data sequence. It is a measure of the average size of "
				+ "the elements of a data sequence.",
		comment = "",
		examples = { @example (" list<float> data_sequence <- [6.0, 7.0, 8.0, 9.0]; "),
				@example (" list<float> squares <- data_sequence collect (each*each); "), @example (
						value = " rms(length(data_sequence),sum(squares)) with_precision(4) ",
						equals = "7.5829") })
public static Double opRms(final IScope scope, final Integer size, final Double sumOfSquares) {

	// TODO input parameters validation

	return Descriptive.rms(size, sumOfSquares);
}
 

/**
 *
 *
 * @param scope
 * @param data
 * @param element
 * @return
 */
@operator (
		value = "rank_interpolated",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the linearly interpolated number of elements in a list less or equal to a given element. The rank is the number of elements <= element. Ranks are of the form {0, 1, 2,..., sortedList.size()}. If no element is <= element, then the rank is zero. If the element lies in between two contained elements, then linear interpolation is used and a non integer value is returned. Note that the container holding the values must be sorted first",
		comment = "",
		examples = { @example (
				value = "rank_interpolated([1,3,5,6,9,11,12,13,19,21,22,32,35,36,45,44,55,68,79,80,81,88,90,91,92,100], 35)",
				equals = "13.0") })
public static Double opRankInterpolated(final IScope scope, final IContainer data, final Double element) {

	// TODO input parameters validation

	return Descriptive.rankInterpolated(from(scope, data), element);
}
 

/**
 *
 *
 * @param scope
 * @param data
 * @param phi
 * @return
 */
@operator (
		value = "quantile",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the phi-quantile; that is, an element elem for which holds that phi percent of data elements are less than elem. The quantile need not necessarily be contained in the data sequence, it can be a linear interpolation. Note that the container holding the values must be sorted first",
		comment = "",
		examples = { @example (
				value = "quantile([1,3,5,6,9,11,12,13,19,21,22,32,35,36,45,44,55,68,79,80,81,88,90,91,92,100], 0.5)",
				equals = "35.5") })
public static Double opQuantile(final IScope scope, final IContainer data, final Double phi) {

	// TODO input parameters validation

	return Descriptive.quantile(from(scope, data), phi);
}
 

/**
 * Returns the moment of k-th order with constant c of a data sequence, which is Sum( (data[i]-c)k ) /
 * data.size().
 *
 * @param scope
 * @param data
 * @param k
 * @param c
 * @return
 */
@operator (
		value = "moment",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the moment of k-th order with constant c of a data sequence",
		comment = "",
		examples = { @example (
				value = "moment([13,2,1,4,1,2], 2, 1.2) with_precision(4)",
				equals = "24.74") })
public static Double opMoment(final IScope scope, final IContainer data, final Integer k, final Double c) {

	// TODO input parameters validation

	return Descriptive.moment(from(scope, data), k, c);
}
 

/**
 *
 *
 * @param scope
 * @param moment4
 * @param standardDeviation
 * @return
 */
@operator (
		value = "kurtosis",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the kurtosis (aka excess) of a data sequence",
		comment = "",
		examples = { @example (
				value = "kurtosis(3,12) with_precision(4)",
				equals = "-2.9999") })
public static Double opKurtosis(final IScope scope, final Double moment4, final Double standardDeviation) {

	// TODO input parameters validation

	return Descriptive.kurtosis(moment4, standardDeviation);
}
 

/**
 *
 *
 * @param scope
 * @param data
 * @return
 */
@operator (
		value = "kurtosis",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the kurtosis (aka excess) of a data sequence",
		comment = "",
		examples = { @example (
				value = "kurtosis([13,2,1,4,1,2]) with_precision(4)",
				equals = "4.8083") })
public static Double opKurtosis(final IScope scope, final IContainer data) {

	// TODO input parameters validation

	final double mean = (Double) Containers.mean(scope, data);
	final double standardDeviation = Stats.opStDev(scope, data);

	return Descriptive.kurtosis(from(scope, data), mean, standardDeviation);
}
 

/**
 *
 *
 * @param scope
 * @param data
 * @return
 */
@operator (
		value = "durbin_watson",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Durbin-Watson computation",
		comment = "",
		examples = { @example (
				value = "durbin_watson([13,2,1,4,1,2]) with_precision(4)",
				equals = "0.7231") })
public static Double opDurbinWatson(final IScope scope, final IContainer data) {

	// TODO input parameters validation

	return Descriptive.durbinWatson(from(scope, data));
}
 

/**
 *
 *
 * @param scope
 * @param data1
 * @param data2
 * @return
 */
@operator (
		value = "covariance",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the covariance of two data sequences",
		comment = "",
		examples = { @example (
				value = "covariance([13,2,1,4,1,2], [1,2,1,3,1,2]) with_precision(2)",
				equals = "-0.67") })
public static Double opCovariance(final IScope scope, final IContainer data1, final IContainer data2) {

	// TODO input parameters validation

	return Descriptive.covariance(from(scope, data1), from(scope, data2));
}
 

/**
 * Returns the correlation of two data sequences.
 *
 * @see <a href="http://www.mathsisfun.com/data/correlation.html"> Correlation</a>
 *
 * @param scope
 * @param data1
 * @param standardDev1
 * @param data2
 * @param stanardDev2
 * @return
 */
@operator (
		value = "correlation",
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns the correlation of two data sequences (having the same size)",
		comment = "",
		examples = { @example (
				value = "correlation([1,2,1,3,1,2], [1,2,1,3,1,2]) with_precision(4)",
				equals = "1.2"),
				@example (
						value = "correlation([13,2,1,4,1,2], [1,2,1,3,1,2]) with_precision(2)",
						equals = "-0.21") })
public static Double opCorrelation(final IScope scope, final IContainer data1, final IContainer data2) {

	// TODO input parameters validation

	final double standardDev1 = Stats.opStDev(scope, data1);
	final double standardDev2 = Stats.opStDev(scope, data2);

	return Descriptive.correlation(from(scope, data1), standardDev1, from(scope, data2), standardDev2);
}
 

/**
 *
 *
 * @param scope
 * @param data
 * @param element
 * @return
 */
@operator (
		value = { "quantile_inverse", "percentile" },
		can_be_const = true,
		type = IType.FLOAT,
		expected_content_type = { IType.INT, IType.FLOAT },
		category = { IOperatorCategory.STATISTICAL },
		concept = { IConcept.STATISTIC })
@doc (
		value = "Returns how many percent of the elements contained in the receiver are <= element. Does linear interpolation if the element is not contained but lies in between two contained elements. Note that the container holding the values must be sorted first",
		comment = "",
		examples = { @example (
				value = "quantile_inverse([1,3,5,6,9,11,12,13,19,21,22,32,35,36,45,44,55,68,79,80,81,88,90,91,92,100], 35.5) with_precision(2)",
				equals = "0.52") })
public static Double opQuantileInverse(final IScope scope, final IContainer data, final Double element) {

	// TODO input parameters validation

	return Descriptive.quantileInverse(from(scope, data), element);
}
 

/**
 * Modifies the receiver to be standardized.
 * Changes each element <tt>x[i]</tt> as follows: <tt>x[i] = (x[i]-mean)/standardDeviation</tt>.
 */
public synchronized void standardize(double mean, double standardDeviation) {
	Descriptive.standardize(this.elements, mean, standardDeviation);
	clearAllMeasures();
	invalidateAll();
	this.size = 0;
}
 

/**
 * Returns the <tt>k-th</tt> order sum of powers, which is <tt>Sum( x[i]<sup>k</sup> )</tt>.
 * @param k the order of the powers.
 * @return the sum of powers.
 */
public synchronized double sumOfPowers(int k) {
	// no chaching for this measure
	if (k >= -1 && k <= 2) return super.sumOfPowers(k);

	return Descriptive.sumOfPowers(this.elements,k);
}
 

/**
 * Returns the moment of <tt>k</tt>-th order with value <tt>c</tt>,
 * which is <tt>Sum( (x[i]-c)<sup>k</sup> ) / size()</tt>.
 *
 * @param k the order; must be greater than or equal to zero.
 * @param c any number.
 * @throws IllegalArgumentException if <tt>k < 0</tt>.
 * @return <tt>Double.NaN</tt> if <tt>!hasSumOfPower(k)</tt>.
 */
public synchronized double moment(int k, double c) {
	if (k<0) throw new IllegalArgumentException("k must be >= 0");
	//checkOrder(k);
	if (!hasSumOfPowers(k)) return Double.NaN;

	int maxOrder = Math.min(k,getMaxOrderForSumOfPowers());
	DoubleArrayList sumOfPows = new DoubleArrayList(maxOrder+1);
	sumOfPows.add(size());
	sumOfPows.add(sum());
	sumOfPows.add(sumOfSquares());
	for (int i=3; i<=maxOrder; i++) sumOfPows.add(sumOfPowers(i));
	
	return Descriptive.moment(k, c, size(), sumOfPows.elements());
}
 

public double bandwidthNRD(double[] in) {

        DoubleArrayList inList = new DoubleArrayList(in.length);
        for (double d : in)
            inList.add(d);
        inList.sort();

        final double h = (Descriptive.quantile(inList, 0.75) - Descriptive.quantile(inList, 0.25)) / 1.34;

        return 4 * 1.06 *
                Math.min(Math.sqrt(DiscreteStatistics.variance(in)), h) *
                Math.pow(in.length, -0.2);
    }
 

/**
 * Returns the rms (Root Mean Square), which is <tt>Math.sqrt( Sum( x[i]*x[i] ) / size() )</tt>.
 */
public synchronized double rms() {
	return Descriptive.rms(size(), sumOfSquares());
}
 

/**
 * 
 */
protected void updateSumOfLogarithms() {
	this.sumOfLogarithms = Descriptive.sumOfLogarithms(this.elements,0,size()-1);
	this.isSumOfLogarithmsValid = true;
}
 

/**
 * assertion: isBasicParametersValid == false
 * 
 */
protected void updateSumOfInversions() {
	this.sumOfInversions = Descriptive.sumOfInversions(this.elements,0,size()-1);
	this.isSumOfInversionsValid = true;
}
 

/**
 * Returns the sample variance, which is <tt>Sum( (x[i]-mean())<sup>2</sup> )  /  (size()-1)</tt>.
 */
public synchronized double variance() {
	return Descriptive.sampleVariance(size(), sum(), sumOfSquares());
}
 

/**
 * Returns the sample standard error, which is <tt>Math.sqrt(variance() / size())</tt>
 */
public synchronized double standardError() {
	return Descriptive.standardError(size(),variance());
}
 

/**
 * Returns the exact quantiles of the specified percentages.
 * @param percentages the percentages for which quantiles are to be computed.
 * Each percentage must be in the interval <tt>(0.0,1.0]</tt>. <tt>percentages</tt> must be sorted ascending.
 * @return the exact quantiles.
 */
public DoubleArrayList quantiles(DoubleArrayList percentages) {
	return Descriptive.quantiles(sortedElements_unsafe(),percentages);
}
 

/**
 * Returns exactly how many percent of the elements contained in the receiver are <tt>&lt;= element</tt>.
 * Does linear interpolation if the element is not contained but lies in between two contained elements.
 *
 * @param element the element to search for.
 * @return the exact percentage <tt>phi</tt> of elements <tt>&lt;= element</tt> (<tt>0.0 &lt;= phi &lt;= 1.0)</tt>.
 */
public synchronized double quantileInverse(double element) {
	return Descriptive.quantileInverse(sortedElements_unsafe(),element);
}
 

/**
 * Returns the exact <tt>phi-</tt>quantile; that is, the smallest contained element <tt>elem</tt> for which holds that <tt>phi</tt> percent of elements are less than <tt>elem</tt>.
 * @param phi must satisfy <tt>0 &lt; phi &lt; 1</tt>.
 */
public synchronized double quantile(double phi) {
	return Descriptive.quantile(sortedElements_unsafe(),phi);
}