Python TsFresh特征的Java实现——max_langevin_fixed_point

TsFresh（TimeSeries Fresh）是一个Python第三方工具包。它可以方便地对时间序列数据进行处理，获得大量的特征。这些特征可以用以训练分类器，以高效地实现对时间序列数据的分类、识别等。然而，在工程实现时，更多地是采用Java等语言，这需要利用Java实现对TsFresh的特征进行直接计算，故需要对TsFresh的某些特征进行深入地分析，并在Java语言下实现。

特征max_langevin_fixed_point简介

命令格式：max_langevin_fixed_point(x, r, m)

根据表征Langevin模型的确定动力学的多项式h(x)估计出动力学的最大固定点

参数：

• x：时间序列，数据类型：numpy.ndarray
• r：分位数的数量，数据类型：float
• m：拟合多项式的除数，数据类型：int

特征max_langevin_fixed_point计算原理

     1. 构造由原时间序列和时间序列x的一阶差分div构成的数据矩阵x_data，第一列为原时间序列，第二列为一阶差分。由于一阶差分div的数据长度较时间序列x的长度少1，故数据矩阵x_data的行数按div的长度而定

     2. 对数据矩阵x_data按第一列的数据由小到大进行排序

Arrays.sort(x_data, new Comparator<double[]>() {
    @Override
    public int compare(double[] t0, double[] t1) { // t0与t1类似于两个对象进行比较，这里表示两个一维数组，即二维数组的两行
        return (int) (t0[0] - t1[0]);
    }
});

     3. 根据数据矩阵x_data第一列进行分位数计算，形成分位数构成的统计区间

     4. 根据数据矩阵x_data第一列的数据所在的分位数区间，为数据矩阵x_data每一行数据分配相应的区间

     5. 统计各个分配数据的统计区间的数据均值，及相应的一阶差分的均值，并赋值给WeightedObservedPoints obs变量

     6. 定义m阶拟合多项式：

PolynomialCurveFitter fitter = PolynomialCurveFitter.create(m);

     7. 根据obs求解m阶拟合多项式的系数：

double[] coeff = fitter.fit(obs.toList());

     8. 求解多项式的根

LaguerreSolver laguerreSolver = new LaguerreSolver();
Complex[] roots = laguerreSolver.solveAllComplex(coeff, -1000000);

     9. 选出最大实部

Java实现

需要用到下列类

import org.apache.commons.math3.analysis.solvers.LaguerreSolver;
import org.apache.commons.math3.fitting.PolynomialCurveFitter;
import org.apache.commons.math3.fitting.WeightedObservedPoints;
import org.apache.commons.math3.complex.Complex;

实现max_langevin_fixed_point特征的代码

	public static double max_langevin_fixed_point(double[] x, double r, int m) {
		int n = x.length;
		double[][] x_data = new double[n - 1][3];
		for (int idx = 0; idx < (n-1); idx++) {
			x_data[idx][0] = x[idx];
			x_data[idx][1] = x[idx + 1] - x[idx];
			//System.out.println("[" + x_data[idx][0] + "," + x_data[idx][1] + "," + x_data[idx][2] + "]");
		}
		//System.out.println("===========================================");
	     //x_data以第一列值大小排序
		Arrays.sort(x_data, new Comparator<double[]>() {
           @Override
           public int compare(double[] t0, double[] t1) { // t0与t1类似于两个对象进行比较，这里表示两个一维数组，即二维数组的两行
               //if(t0[1] == t1[1]){
               //    return (int) (t0[0]-t1[0]); //第二列中的值相同时，再比较第一列的大小，如果不需要也可以直接return 0
               //}
               return (int) (t0[0] - t1[0]);
           }
       });
		for (int idx = 0; idx < (n-1); idx++) {
			System.out.println("[" + x_data[idx][0] + "," + x_data[idx][1] + "," + x_data[idx][2] + "]");
		}
		
		// 根据x_data的第一列数据，按r分位数将数据区间进行分割区间
		double[] x_sorted = new double[n-1];
		for (int idx = 0; idx < (n-1); idx++) {
			x_sorted[idx] = x_data[idx][0];
		}
		// blocks用于表示分位数的值，其与x_sorted[0]一起将数据分割为r段
		double[] blocks = new double[(int) r];
		for (int idx = 0; idx < blocks.length; idx++) {
			blocks[idx] = Quantiles_INC(x_sorted, ((double) (idx +1))/r);
			//System.out.println(blocks[idx]);
		}
		// 判断x_data所处的区间
		int index = 0;
		x_data[0][2] = blocks[index];
		for (int idx = 1; idx < (n-1); idx++) {
			for (int ind = index; ind < blocks.length; ind++) {
				if (ind == 0) {
					if (x_data[idx][0] <= blocks[ind]) {
						x_data[idx][2] = blocks[ind];
						break;
					}
				} else {
					if ((x_data[idx][0] <= blocks[ind])&(x_data[idx][0] > blocks[ind - 1])) {
						x_data[idx][2] = blocks[ind];
						index = ind;
						break;
					}
				}
			}
		}
		for (int idx = 0; idx < (n-1); idx++) {
			System.out.println("[" + x_data[idx][0] + "," + x_data[idx][1] + "," + x_data[idx][2] + "]");
		}
		
		// 求每个区间的均值
		double tempData0 = x_data[0][0];
		double tempData1 = x_data[0][1];
		int tempNum = 1;
		index = 0;
		WeightedObservedPoints obs = new WeightedObservedPoints();
		for (int idx = 1; idx <(n-1); idx++) {
			if (x_data[idx][2] > x_data[idx - 1][2]) {
				obs.add(tempData0/tempNum, tempData1/tempNum);
				//System.out.println("[" + tempData0/tempNum + "\t" + tempData1/tempNum + "]");
				tempData0 = x_data[idx][0];
				tempData1 = x_data[idx][1];
				tempNum = 1;
				index++;
			} else {
				tempData0 = tempData0 + x_data[idx][0];
				tempData1 = tempData1 + x_data[idx][1];
				tempNum++;
			}
			if (idx == (n-2)) {
				obs.add(tempData0/tempNum, tempData1/tempNum);
				//System.out.println("[" + tempData0/tempNum + "\t" + tempData1/tempNum + "]");
			}

		}
		// Instantiate a m-degree polynomial fitter.
		final PolynomialCurveFitter fitter = PolynomialCurveFitter.create(m);

		// Retrieve fitted parameters (coefficients of the polynomial function).
		final double[] coeff = fitter.fit(obs.toList());
		//System.out.println(Arrays.toString(coeff));
		
		// 求根
		LaguerreSolver laguerreSolver = new LaguerreSolver();
		Complex[] xx = laguerreSolver.solveAllComplex(coeff, -1000000);
		//System.out.println("全部复根：" + xx[0] + "\t" + xx[1] + "\t" + xx[2]);
		// 选择出最大实部
		double result = xx[0].getReal();
		for (int idx = 1; idx < m; idx++) {
			if (xx[idx].getReal() > result) {
				result = xx[idx].getReal(); 
			}
		}
		
		return result;
	}