source: branches/3040_VectorBasedGP/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/SymbolicDataAnalysisExpressionTreeVectorInterpreter.cs @ 18058

Last change on this file since 18058 was 18058, checked in by pfleck, 11 months ago

#3040 Switched from offset-length sub-vector to "start-end with wrapping" subvector.

File size: 51.4 KB
Line 
1#region License Information
2/* HeuristicLab
3 * Copyright (C) Heuristic and Evolutionary Algorithms Laboratory (HEAL)
4 *
5 * This file is part of HeuristicLab.
6 *
7 * HeuristicLab is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * HeuristicLab is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
19 */
20#endregion
21
22using System;
23using System.Collections.Generic;
24using System.Linq;
25using HeuristicLab.Common;
26using HeuristicLab.Core;
27using HeuristicLab.Data;
28using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
29using HeuristicLab.Parameters;
30using HEAL.Attic;
31using MathNet.Numerics;
32using MathNet.Numerics.Statistics;
33using DoubleVector = MathNet.Numerics.LinearAlgebra.Vector<double>;
34
35namespace HeuristicLab.Problems.DataAnalysis.Symbolic {
36  [StorableType("DE68A1D9-5AFC-4DDD-AB62-29F3B8FC28E0")]
37  [Item("SymbolicDataAnalysisExpressionTreeVectorInterpreter", "Interpreter for symbolic expression trees including vector arithmetic.")]
38  public class SymbolicDataAnalysisExpressionTreeVectorInterpreter : ParameterizedNamedItem, ISymbolicDataAnalysisExpressionTreeInterpreter {
39    [StorableType("2612504E-AD5F-4AE2-B60E-98A5AB59E164")]
40    public enum Aggregation {
41      Mean,
42      Median,
43      Sum,
44      First,
45      L1Norm,
46      L2Norm,
47      NaN,
48      Exception
49    }
50    public static double Aggregate(Aggregation aggregation, DoubleVector vector) {
51      switch (aggregation) {
52        case Aggregation.Mean: return Statistics.Mean(vector);
53        case Aggregation.Median: return Statistics.Median(vector);
54        case Aggregation.Sum: return vector.Sum();
55        case Aggregation.First: return vector.First();
56        case Aggregation.L1Norm: return vector.L1Norm();
57        case Aggregation.L2Norm: return vector.L2Norm();
58        case Aggregation.NaN: return double.NaN;
59        case Aggregation.Exception: throw new InvalidOperationException("Result of the tree is not a scalar.");
60        default: throw new ArgumentOutOfRangeException(nameof(aggregation), aggregation, null);
61      }
62    }
63
64    [StorableType("73DCBB45-916F-4139-8ADC-57BA610A1B66")]
65    public enum VectorLengthStrategy {
66      ExceptionIfDifferent,
67      FillShorterWithNaN,
68      FillShorterWithNeutralElement,
69      CutLonger,
70      ResampleToLonger,
71      ResampleToShorter,
72      CycleShorter
73    }
74
75    #region Implementation VectorLengthStrategy
76    public static (DoubleVector, DoubleVector) ExceptionIfDifferent(DoubleVector lhs, DoubleVector rhs) {
77      if (lhs.Count != rhs.Count)
78        throw new InvalidOperationException($"Vector Lengths incompatible ({lhs.Count} vs. {rhs.Count}");
79      return (lhs, rhs);
80    }
81
82    public static (DoubleVector, DoubleVector) FillShorter(DoubleVector lhs, DoubleVector rhs, double fillElement) {
83      var targetLength = Math.Max(lhs.Count, rhs.Count);
84
85      DoubleVector PadVector(DoubleVector v) {
86        if (v.Count == targetLength) return v;
87        var p = DoubleVector.Build.Dense(targetLength, fillElement);
88        v.CopySubVectorTo(p, 0, 0, v.Count);
89        return p;
90      }
91
92      return (PadVector(lhs), PadVector(rhs));
93    }
94
95    public static (DoubleVector, DoubleVector) CutLonger(DoubleVector lhs, DoubleVector rhs) {
96      var targetLength = Math.Min(lhs.Count, rhs.Count);
97
98      DoubleVector CutVector(DoubleVector v) {
99        if (v.Count == targetLength) return v;
100        return v.SubVector(0, targetLength);
101      }
102
103      return (CutVector(lhs), CutVector(rhs));
104    }
105
106    private static DoubleVector ResampleToLength(DoubleVector v, int targetLength) {
107      if (v.Count == targetLength) return v;
108
109      var indices = Enumerable.Range(0, v.Count).Select(x => (double)x);
110      var interpolation = Interpolate.Linear(indices, v);
111
112      var resampledIndices = Enumerable.Range(0, targetLength).Select(i => (double)i / targetLength * v.Count);
113      var interpolatedValues = resampledIndices.Select(interpolation.Interpolate);
114
115      return DoubleVector.Build.DenseOfEnumerable(interpolatedValues);
116    }
117    public static (DoubleVector, DoubleVector) ResampleToLonger(DoubleVector lhs, DoubleVector rhs) {
118      var maxLength = Math.Max(lhs.Count, rhs.Count);
119      return (ResampleToLength(lhs, maxLength), ResampleToLength(rhs, maxLength));
120    }
121    public static (DoubleVector, DoubleVector) ResampleToShorter(DoubleVector lhs, DoubleVector rhs) {
122      var minLength = Math.Min(lhs.Count, rhs.Count);
123      return (ResampleToLength(lhs, minLength), ResampleToLength(rhs, minLength));
124    }
125
126    public static (DoubleVector, DoubleVector) CycleShorter(DoubleVector lhs, DoubleVector rhs) {
127      var targetLength = Math.Max(lhs.Count, rhs.Count);
128
129      DoubleVector CycleVector(DoubleVector v) {
130        if (v.Count == targetLength) return v;
131        var cycledValues = Enumerable.Range(0, targetLength).Select(i => v[i % v.Count]);
132        return DoubleVector.Build.DenseOfEnumerable(cycledValues);
133      }
134
135      return (CycleVector(lhs), CycleVector(rhs));
136    }
137    #endregion
138
139    public static (DoubleVector lhs, DoubleVector rhs) ApplyVectorLengthStrategy(VectorLengthStrategy strategy, DoubleVector lhs, DoubleVector rhs,
140      double neutralElement = double.NaN) {
141
142      switch (strategy) {
143        case VectorLengthStrategy.ExceptionIfDifferent: return ExceptionIfDifferent(lhs, rhs);
144        case VectorLengthStrategy.FillShorterWithNaN: return FillShorter(lhs, rhs, double.NaN);
145        case VectorLengthStrategy.FillShorterWithNeutralElement: return FillShorter(lhs, rhs, neutralElement);
146        case VectorLengthStrategy.CutLonger: return CutLonger(lhs, rhs);
147        case VectorLengthStrategy.ResampleToLonger: return ResampleToLonger(lhs, rhs);
148        case VectorLengthStrategy.ResampleToShorter: return ResampleToShorter(lhs, rhs);
149        case VectorLengthStrategy.CycleShorter: return CycleShorter(lhs, rhs);
150        default: throw new ArgumentOutOfRangeException(nameof(strategy), strategy, null);
151      }
152    }
153
154    #region Aggregation Symbols
155    private static Type[] AggregationSymbols = new[] {
156      typeof(Sum), typeof(Mean), typeof(Length), typeof(StandardDeviation), typeof(Variance),
157      typeof(EuclideanDistance), typeof(Covariance)
158    };
159    #endregion
160
161    private const string EvaluatedSolutionsParameterName = "EvaluatedSolutions";
162    private const string FinalAggregationParameterName = "FinalAggregation";
163    private const string DifferentVectorLengthStrategyParameterName = "DifferentVectorLengthStrategy";
164
165    public override bool CanChangeName {
166      get { return false; }
167    }
168
169    public override bool CanChangeDescription {
170      get { return false; }
171    }
172
173    #region parameter properties
174    public IFixedValueParameter<IntValue> EvaluatedSolutionsParameter {
175      get { return (IFixedValueParameter<IntValue>)Parameters[EvaluatedSolutionsParameterName]; }
176    }
177    public IFixedValueParameter<EnumValue<Aggregation>> FinalAggregationParameter {
178      get { return (IFixedValueParameter<EnumValue<Aggregation>>)Parameters[FinalAggregationParameterName]; }
179    }
180    public IFixedValueParameter<EnumValue<VectorLengthStrategy>> DifferentVectorLengthStrategyParameter {
181      get { return (IFixedValueParameter<EnumValue<VectorLengthStrategy>>)Parameters[DifferentVectorLengthStrategyParameterName]; }
182    }
183    #endregion
184
185    #region properties
186    public int EvaluatedSolutions {
187      get { return EvaluatedSolutionsParameter.Value.Value; }
188      set { EvaluatedSolutionsParameter.Value.Value = value; }
189    }
190    public Aggregation FinalAggregation {
191      get { return FinalAggregationParameter.Value.Value; }
192      set { FinalAggregationParameter.Value.Value = value; }
193    }
194    public VectorLengthStrategy DifferentVectorLengthStrategy {
195      get { return DifferentVectorLengthStrategyParameter.Value.Value; }
196      set { DifferentVectorLengthStrategyParameter.Value.Value = value; }
197    }
198    #endregion
199
200    [StorableConstructor]
201    protected SymbolicDataAnalysisExpressionTreeVectorInterpreter(StorableConstructorFlag _) : base(_) { }
202
203    protected SymbolicDataAnalysisExpressionTreeVectorInterpreter(SymbolicDataAnalysisExpressionTreeVectorInterpreter original, Cloner cloner)
204      : base(original, cloner) { }
205
206    public override IDeepCloneable Clone(Cloner cloner) {
207      return new SymbolicDataAnalysisExpressionTreeVectorInterpreter(this, cloner);
208    }
209
210    public SymbolicDataAnalysisExpressionTreeVectorInterpreter()
211      : this("SymbolicDataAnalysisExpressionTreeVectorInterpreter", "Interpreter for symbolic expression trees including vector arithmetic.") { }
212
213    protected SymbolicDataAnalysisExpressionTreeVectorInterpreter(string name, string description)
214      : base(name, description) {
215      Parameters.Add(new FixedValueParameter<IntValue>(EvaluatedSolutionsParameterName, "A counter for the total number of solutions the interpreter has evaluated", new IntValue(0)));
216      Parameters.Add(new FixedValueParameter<EnumValue<Aggregation>>(FinalAggregationParameterName, "If root node of the expression tree results in a Vector it is aggregated according to this parameter", new EnumValue<Aggregation>(Aggregation.Mean)));
217      Parameters.Add(new FixedValueParameter<EnumValue<VectorLengthStrategy>>(DifferentVectorLengthStrategyParameterName, "", new EnumValue<VectorLengthStrategy>(VectorLengthStrategy.ExceptionIfDifferent)));
218    }
219
220    [StorableHook(HookType.AfterDeserialization)]
221    private void AfterDeserialization() {
222      if (!Parameters.ContainsKey(FinalAggregationParameterName)) {
223        Parameters.Add(new FixedValueParameter<EnumValue<Aggregation>>(FinalAggregationParameterName, "If root node of the expression tree results in a Vector it is aggregated according to this parameter", new EnumValue<Aggregation>(Aggregation.Mean)));
224      }
225      if (!Parameters.ContainsKey(DifferentVectorLengthStrategyParameterName)) {
226        Parameters.Add(new FixedValueParameter<EnumValue<VectorLengthStrategy>>(DifferentVectorLengthStrategyParameterName, "", new EnumValue<VectorLengthStrategy>(VectorLengthStrategy.ExceptionIfDifferent)));
227      }
228    }
229
230    #region IStatefulItem
231    public void InitializeState() {
232      EvaluatedSolutions = 0;
233    }
234
235    public void ClearState() { }
236    #endregion
237
238    private readonly object syncRoot = new object();
239    public IEnumerable<double> GetSymbolicExpressionTreeValues(ISymbolicExpressionTree tree, IDataset dataset, IEnumerable<int> rows) {
240      lock (syncRoot) {
241        EvaluatedSolutions++; // increment the evaluated solutions counter
242      }
243      var state = PrepareInterpreterState(tree, dataset);
244
245      foreach (var rowEnum in rows) {
246        int row = rowEnum;
247        var result = Evaluate(dataset, ref row, state);
248        if (result.IsScalar)
249          yield return result.Scalar;
250        else if (result.IsVector) {
251          yield return Aggregate(FinalAggregation, result.Vector);
252        } else
253          yield return double.NaN;
254        state.Reset();
255      }
256    }
257
258    public IEnumerable<Dictionary<ISymbolicExpressionTreeNode, EvaluationResult>> GetIntermediateNodeValues(ISymbolicExpressionTree tree, IDataset dataset, IEnumerable<int> rows) {
259      var state = PrepareInterpreterState(tree, dataset);
260
261      foreach (var rowEnum in rows) {
262        int row = rowEnum;
263        var traceDict = new Dictionary<ISymbolicExpressionTreeNode, EvaluationResult>();
264        var result = Evaluate(dataset, ref row, state, traceDict);
265        traceDict.Add(tree.Root.GetSubtree(0), result); // Add StartSymbol
266        yield return traceDict;
267        state.Reset();
268      }
269    }
270
271    private static InterpreterState PrepareInterpreterState(ISymbolicExpressionTree tree, IDataset dataset) {
272      Instruction[] code = SymbolicExpressionTreeCompiler.Compile(tree, OpCodes.MapSymbolToOpCode);
273      int necessaryArgStackSize = 0;
274      foreach (Instruction instr in code) {
275        if (instr.opCode == OpCodes.Variable) {
276          var variableTreeNode = (VariableTreeNode)instr.dynamicNode;
277          if (dataset.VariableHasType<double>(variableTreeNode.VariableName))
278            instr.data = dataset.GetReadOnlyDoubleValues(variableTreeNode.VariableName);
279          else if (dataset.VariableHasType<DoubleVector>(variableTreeNode.VariableName))
280            instr.data = dataset.GetReadOnlyDoubleVectorValues(variableTreeNode.VariableName);
281          else throw new NotSupportedException($"Type of variable {variableTreeNode.VariableName} is not supported.");
282        } else if (instr.opCode == OpCodes.FactorVariable) {
283          var factorTreeNode = instr.dynamicNode as FactorVariableTreeNode;
284          instr.data = dataset.GetReadOnlyStringValues(factorTreeNode.VariableName);
285        } else if (instr.opCode == OpCodes.BinaryFactorVariable) {
286          var factorTreeNode = instr.dynamicNode as BinaryFactorVariableTreeNode;
287          instr.data = dataset.GetReadOnlyStringValues(factorTreeNode.VariableName);
288        } else if (instr.opCode == OpCodes.LagVariable) {
289          var laggedVariableTreeNode = (LaggedVariableTreeNode)instr.dynamicNode;
290          instr.data = dataset.GetReadOnlyDoubleValues(laggedVariableTreeNode.VariableName);
291        } else if (instr.opCode == OpCodes.VariableCondition) {
292          var variableConditionTreeNode = (VariableConditionTreeNode)instr.dynamicNode;
293          instr.data = dataset.GetReadOnlyDoubleValues(variableConditionTreeNode.VariableName);
294        } else if (instr.opCode == OpCodes.Call) {
295          necessaryArgStackSize += instr.nArguments + 1;
296        }
297      }
298      return new InterpreterState(code, necessaryArgStackSize);
299    }
300
301
302    public struct EvaluationResult {
303      public double Scalar { get; }
304      public bool IsScalar => !double.IsNaN(Scalar);
305
306      public DoubleVector Vector { get; }
307      public bool IsVector => !(Vector.Count == 1 && double.IsNaN(Vector[0]));
308
309      public bool IsNaN => !IsScalar && !IsVector;
310
311      public EvaluationResult(double scalar) {
312        Scalar = scalar;
313        Vector = NaNVector;
314      }
315      public EvaluationResult(DoubleVector vector) {
316        if (vector == null) throw new ArgumentNullException(nameof(vector));
317        Vector = vector;
318        Scalar = double.NaN;
319      }
320
321      public override string ToString() {
322        if (IsScalar) return Scalar.ToString();
323        if (IsVector) return Vector.ToVectorString();
324        return "NaN";
325      }
326
327      private static readonly DoubleVector NaNVector = DoubleVector.Build.Dense(1, double.NaN);
328      public static readonly EvaluationResult NaN = new EvaluationResult(double.NaN);
329    }
330
331    private static EvaluationResult ArithmeticApply(EvaluationResult lhs, EvaluationResult rhs,
332      Func<DoubleVector, DoubleVector, (DoubleVector, DoubleVector)> lengthStrategy,
333      Func<double, double, double> ssFunc = null,
334      Func<double, DoubleVector, DoubleVector> svFunc = null,
335      Func<DoubleVector, double, DoubleVector> vsFunc = null,
336      Func<DoubleVector, DoubleVector, DoubleVector> vvFunc = null) {
337
338      if (lhs.IsScalar && rhs.IsScalar && ssFunc != null) return new EvaluationResult(ssFunc(lhs.Scalar, rhs.Scalar));
339      if (lhs.IsScalar && rhs.IsVector && svFunc != null) return new EvaluationResult(svFunc(lhs.Scalar, rhs.Vector));
340      if (lhs.IsVector && rhs.IsScalar && vsFunc != null) return new EvaluationResult(vsFunc(lhs.Vector, rhs.Scalar));
341      if (lhs.IsVector && rhs.IsVector && vvFunc != null) {
342        if (lhs.Vector.Count == rhs.Vector.Count) {
343          return new EvaluationResult(vvFunc(lhs.Vector, rhs.Vector));
344        } else {
345          var (lhsVector, rhsVector) = lengthStrategy(lhs.Vector, rhs.Vector);
346          return new EvaluationResult(vvFunc(lhsVector, rhsVector));
347        }
348      }
349      return EvaluationResult.NaN;
350    }
351
352    private static EvaluationResult FunctionApply(EvaluationResult val,
353      Func<double, double> sFunc = null,
354      Func<DoubleVector, DoubleVector> vFunc = null) {
355      if (val.IsScalar && sFunc != null) return new EvaluationResult(sFunc(val.Scalar));
356      if (val.IsVector && vFunc != null) return new EvaluationResult(vFunc(val.Vector));
357      return EvaluationResult.NaN;
358    }
359    private static EvaluationResult AggregateApply(EvaluationResult val,
360      Func<double, double> sFunc = null,
361      Func<DoubleVector, double> vFunc = null) {
362      if (val.IsScalar && sFunc != null) return new EvaluationResult(sFunc(val.Scalar));
363      if (val.IsVector && vFunc != null) return new EvaluationResult(vFunc(val.Vector));
364      return EvaluationResult.NaN;
365    }
366
367    private static EvaluationResult WindowedAggregateApply(EvaluationResult val, WindowedSymbolTreeNode node,
368      Func<double, double> sFunc = null,
369      Func<DoubleVector, double> vFunc = null) {
370
371      // Parameters are interpreted as start and end with wrapping
372      var start = node.Offset;
373      var end = node.Length;
374
375      DoubleVector SubVector(DoubleVector v) {
376        int startIdx = (int)Math.Round(start * v.Count);
377        int endIdx = (int)Math.Round(end * v.Count);
378        int size = v.Count;
379        if (startIdx < endIdx) {
380          return v.SubVector(startIdx, count: endIdx - startIdx);
381        } else { // wrap around
382          var resultVector = DoubleVector.Build.Dense(size: size - (startIdx - endIdx));
383          v.CopySubVectorTo(resultVector, startIdx, 0, size - startIdx); // copy [startIdx:size] to [0:size-startIdx]
384          v.CopySubVectorTo(resultVector, 0, size - startIdx, endIdx);   // copy [0:endIdx]      to [size-startIdx:size]
385          return resultVector;
386        }
387      }
388
389      if (val.IsScalar && sFunc != null) return new EvaluationResult(sFunc(val.Scalar));
390      if (val.IsVector && vFunc != null) return new EvaluationResult(vFunc(SubVector(val.Vector)));
391      return EvaluationResult.NaN;
392    }
393    private static EvaluationResult WindowedFunctionApply(EvaluationResult val, IWindowedSymbolTreeNode node,
394      Func<double, double> sFunc = null,
395      Func<DoubleVector, DoubleVector> vFunc = null) {
396      // Parameters are interpreted as start and end with wrapping
397      var start = node.Offset;
398      var end = node.Length;
399
400      DoubleVector SubVector(DoubleVector v) {
401        int startIdx = (int)Math.Round(start * v.Count);
402        int endIdx = (int)Math.Round(end * v.Count);
403        int size = v.Count;
404        if (startIdx < endIdx) {
405          return v.SubVector(startIdx, count: endIdx - startIdx);
406        } else { // wrap around
407          var resultVector = DoubleVector.Build.Dense(size: size - (startIdx - endIdx));
408          v.CopySubVectorTo(resultVector, startIdx, 0, size - startIdx); // copy [startIdx:size] to [0:size-startIdx]
409          v.CopySubVectorTo(resultVector, 0, size - startIdx, endIdx);   // copy [0:endIdx]      to [size-startIdx:size]
410          return resultVector;
411        }
412      }
413
414      if (val.IsScalar && sFunc != null) return new EvaluationResult(sFunc(val.Scalar));
415      if (val.IsVector && vFunc != null) return new EvaluationResult(vFunc(SubVector(val.Vector)));
416      return EvaluationResult.NaN;
417    }
418
419    private static EvaluationResult AggregateMultipleApply(EvaluationResult lhs, EvaluationResult rhs,
420      Func<DoubleVector, DoubleVector, (DoubleVector, DoubleVector)> lengthStrategy,
421      Func<double, double, double> ssFunc = null,
422      Func<double, DoubleVector, double> svFunc = null,
423      Func<DoubleVector, double, double> vsFunc = null,
424      Func<DoubleVector, DoubleVector, double> vvFunc = null) {
425      if (lhs.IsScalar && rhs.IsScalar && ssFunc != null) return new EvaluationResult(ssFunc(lhs.Scalar, rhs.Scalar));
426      if (lhs.IsScalar && rhs.IsVector && svFunc != null) return new EvaluationResult(svFunc(lhs.Scalar, rhs.Vector));
427      if (lhs.IsVector && rhs.IsScalar && vsFunc != null) return new EvaluationResult(vsFunc(lhs.Vector, rhs.Scalar));
428      if (lhs.IsVector && rhs.IsVector && vvFunc != null) {
429        if (lhs.Vector.Count == rhs.Vector.Count) {
430          return new EvaluationResult(vvFunc(lhs.Vector, rhs.Vector));
431        } else {
432          var (lhsVector, rhsVector) = lengthStrategy(lhs.Vector, rhs.Vector);
433          return new EvaluationResult(vvFunc(lhsVector, rhsVector));
434        }
435      }
436      return EvaluationResult.NaN;
437    }
438
439    public virtual Type GetNodeType(ISymbolicExpressionTreeNode node) {
440      if (node.DataType != null)
441        return node.DataType;
442
443      if (AggregationSymbols.Contains(node.Symbol.GetType()))
444        return typeof(double);
445
446      var argumentTypes = node.Subtrees.Select(GetNodeType);
447      if (argumentTypes.Any(t => t == typeof(DoubleVector)))
448        return typeof(DoubleVector);
449
450      return typeof(double);
451    }
452
453
454    public virtual EvaluationResult Evaluate(IDataset dataset, ref int row, InterpreterState state,
455      IDictionary<ISymbolicExpressionTreeNode, EvaluationResult> traceDict = null) {
456
457      void TraceEvaluation(Instruction instr, EvaluationResult result) {
458        traceDict?.Add(instr.dynamicNode, result);
459      }
460
461      Instruction currentInstr = state.NextInstruction();
462      switch (currentInstr.opCode) {
463        case OpCodes.Add: {
464            var cur = Evaluate(dataset, ref row, state, traceDict);
465            for (int i = 1; i < currentInstr.nArguments; i++) {
466              var op = Evaluate(dataset, ref row, state, traceDict);
467              cur = ArithmeticApply(cur, op,
468                (lhs, rhs) => ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 0.0),
469                (s1, s2) => s1 + s2,
470                (s1, v2) => s1 + v2,
471                (v1, s2) => v1 + s2,
472                (v1, v2) => v1 + v2);
473            }
474            TraceEvaluation(currentInstr, cur);
475            return cur;
476          }
477        case OpCodes.Sub: {
478            var cur = Evaluate(dataset, ref row, state, traceDict);
479            for (int i = 1; i < currentInstr.nArguments; i++) {
480              var op = Evaluate(dataset, ref row, state, traceDict);
481              cur = ArithmeticApply(cur, op,
482                (lhs, rhs) => ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 0.0),
483                (s1, s2) => s1 - s2,
484                (s1, v2) => s1 - v2,
485                (v1, s2) => v1 - s2,
486                (v1, v2) => v1 - v2);
487            }
488            if (currentInstr.nArguments == 1)
489              cur = FunctionApply(cur,
490                s => -s,
491                v => -v);
492            TraceEvaluation(currentInstr, cur);
493            return cur;
494          }
495        case OpCodes.Mul: {
496            var cur = Evaluate(dataset, ref row, state, traceDict);
497            for (int i = 1; i < currentInstr.nArguments; i++) {
498              var op = Evaluate(dataset, ref row, state, traceDict);
499              cur = ArithmeticApply(cur, op,
500                (lhs, rhs) => ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 1.0),
501                (s1, s2) => s1 * s2,
502                (s1, v2) => s1 * v2,
503                (v1, s2) => v1 * s2,
504                (v1, v2) => v1.PointwiseMultiply(v2));
505            }
506            TraceEvaluation(currentInstr, cur);
507            return cur;
508          }
509        case OpCodes.Div: {
510            var cur = Evaluate(dataset, ref row, state, traceDict);
511            for (int i = 1; i < currentInstr.nArguments; i++) {
512              var op = Evaluate(dataset, ref row, state, traceDict);
513              cur = ArithmeticApply(cur, op,
514                (lhs, rhs) => ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 1.0),
515                (s1, s2) => s1 / s2,
516                (s1, v2) => s1 / v2,
517                (v1, s2) => v1 / s2,
518                (v1, v2) => v1 / v2);
519            }
520            if (currentInstr.nArguments == 1)
521              cur = FunctionApply(cur,
522                s => 1 / s,
523                v => 1 / v);
524            TraceEvaluation(currentInstr, cur);
525            return cur;
526          }
527        case OpCodes.Absolute: {
528            var cur = Evaluate(dataset, ref row, state, traceDict);
529            cur = FunctionApply(cur, Math.Abs, DoubleVector.Abs);
530            TraceEvaluation(currentInstr, cur);
531            return cur;
532          }
533        case OpCodes.Tanh: {
534            var cur = Evaluate(dataset, ref row, state, traceDict);
535            cur = FunctionApply(cur, Math.Tanh, DoubleVector.Tanh);
536            TraceEvaluation(currentInstr, cur);
537            return cur;
538          }
539        case OpCodes.Cos: {
540            var cur = Evaluate(dataset, ref row, state, traceDict);
541            cur = FunctionApply(cur, Math.Cos, DoubleVector.Cos);
542            TraceEvaluation(currentInstr, cur);
543            return cur;
544          }
545        case OpCodes.Sin: {
546            var cur = Evaluate(dataset, ref row, state, traceDict);
547            cur = FunctionApply(cur, Math.Sin, DoubleVector.Sin);
548            TraceEvaluation(currentInstr, cur);
549            return cur;
550          }
551        case OpCodes.Tan: {
552            var cur = Evaluate(dataset, ref row, state, traceDict);
553            cur = FunctionApply(cur, Math.Tan, DoubleVector.Tan);
554            TraceEvaluation(currentInstr, cur);
555            return cur;
556          }
557        case OpCodes.Square: {
558            var cur = Evaluate(dataset, ref row, state, traceDict);
559            cur = FunctionApply(cur,
560              s => Math.Pow(s, 2),
561              v => v.PointwisePower(2));
562            TraceEvaluation(currentInstr, cur);
563            return cur;
564          }
565        case OpCodes.Cube: {
566            var cur = Evaluate(dataset, ref row, state, traceDict);
567            cur = FunctionApply(cur,
568              s => Math.Pow(s, 3),
569              v => v.PointwisePower(3));
570            TraceEvaluation(currentInstr, cur);
571            return cur;
572          }
573        case OpCodes.Power: {
574            var x = Evaluate(dataset, ref row, state, traceDict);
575            var y = Evaluate(dataset, ref row, state, traceDict);
576            var cur = ArithmeticApply(x, y,
577              (lhs, rhs) => lhs.Count < rhs.Count
578                ? CutLonger(lhs, rhs)
579                : ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 1.0),
580              (s1, s2) => Math.Pow(s1, Math.Round(s2)),
581              (s1, v2) => DoubleVector.Build.Dense(v2.Count, s1).PointwisePower(DoubleVector.Round(v2)),
582              (v1, s2) => v1.PointwisePower(Math.Round(s2)),
583              (v1, v2) => v1.PointwisePower(DoubleVector.Round(v2)));
584            TraceEvaluation(currentInstr, cur);
585            return cur;
586          }
587        case OpCodes.SquareRoot: {
588            var cur = Evaluate(dataset, ref row, state, traceDict);
589            cur = FunctionApply(cur,
590              s => Math.Sqrt(s),
591              v => DoubleVector.Sqrt(v));
592            TraceEvaluation(currentInstr, cur);
593            return cur;
594          }
595        case OpCodes.CubeRoot: {
596            var cur = Evaluate(dataset, ref row, state, traceDict);
597            cur = FunctionApply(cur,
598              s => s < 0 ? -Math.Pow(-s, 1.0 / 3.0) : Math.Pow(s, 1.0 / 3.0),
599              v => v.Map(s => s < 0 ? -Math.Pow(-s, 1.0 / 3.0) : Math.Pow(s, 1.0 / 3.0)));
600            TraceEvaluation(currentInstr, cur);
601            return cur;
602          }
603        case OpCodes.Root: {
604            var x = Evaluate(dataset, ref row, state, traceDict);
605            var y = Evaluate(dataset, ref row, state, traceDict);
606            var cur = ArithmeticApply(x, y,
607              (lhs, rhs) => lhs.Count < rhs.Count
608                ? CutLonger(lhs, rhs)
609                : ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 1.0),
610              (s1, s2) => Math.Pow(s1, 1.0 / Math.Round(s2)),
611              (s1, v2) => DoubleVector.Build.Dense(v2.Count, s1).PointwisePower(1.0 / DoubleVector.Round(v2)),
612              (v1, s2) => v1.PointwisePower(1.0 / Math.Round(s2)),
613              (v1, v2) => v1.PointwisePower(1.0 / DoubleVector.Round(v2)));
614            TraceEvaluation(currentInstr, cur);
615            return cur;
616          }
617        case OpCodes.Exp: {
618            var cur = Evaluate(dataset, ref row, state, traceDict);
619            cur = FunctionApply(cur,
620              s => Math.Exp(s),
621              v => DoubleVector.Exp(v));
622            TraceEvaluation(currentInstr, cur);
623            return cur;
624          }
625        case OpCodes.Log: {
626            var cur = Evaluate(dataset, ref row, state, traceDict);
627            cur = FunctionApply(cur,
628              s => Math.Log(s),
629              v => DoubleVector.Log(v));
630            TraceEvaluation(currentInstr, cur);
631            return cur;
632          }
633        case OpCodes.Sum: {
634            var cur = Evaluate(dataset, ref row, state, traceDict);
635            cur = AggregateApply(cur,
636              s => s,
637              v => v.Sum());
638            TraceEvaluation(currentInstr, cur);
639            return cur;
640          }
641        case OpCodes.Mean: {
642            var cur = Evaluate(dataset, ref row, state, traceDict);
643            cur = AggregateApply(cur,
644              s => s,
645              v => Statistics.Mean(v));
646            TraceEvaluation(currentInstr, cur);
647            return cur;
648          }
649        case OpCodes.StandardDeviation: {
650            var cur = Evaluate(dataset, ref row, state, traceDict);
651            cur = AggregateApply(cur,
652              s => 0,
653              v => Statistics.PopulationStandardDeviation(v));
654            TraceEvaluation(currentInstr, cur);
655            return cur;
656          }
657        case OpCodes.Length: {
658            var cur = Evaluate(dataset, ref row, state, traceDict);
659            cur = AggregateApply(cur,
660              s => 1,
661              v => v.Count);
662            TraceEvaluation(currentInstr, cur);
663            return cur;
664          }
665        case OpCodes.Min: {
666            var cur = Evaluate(dataset, ref row, state, traceDict);
667            cur = AggregateApply(cur,
668              s => s,
669              v => Statistics.Minimum(v));
670            TraceEvaluation(currentInstr, cur);
671            return cur;
672          }
673        case OpCodes.Max: {
674            var cur = Evaluate(dataset, ref row, state, traceDict);
675            cur = AggregateApply(cur,
676              s => s,
677              v => Statistics.Maximum(v));
678            TraceEvaluation(currentInstr, cur);
679            return cur;
680          }
681        case OpCodes.Variance: {
682            var cur = Evaluate(dataset, ref row, state, traceDict);
683            cur = AggregateApply(cur,
684              s => 0,
685              v => Statistics.PopulationVariance(v));
686            TraceEvaluation(currentInstr, cur);
687            return cur;
688          }
689        case OpCodes.Skewness: {
690            var cur = Evaluate(dataset, ref row, state, traceDict);
691            cur = AggregateApply(cur,
692              s => double.NaN,
693              v => Statistics.PopulationSkewness(v));
694            TraceEvaluation(currentInstr, cur);
695            return cur;
696          }
697        case OpCodes.Kurtosis: {
698            var cur = Evaluate(dataset, ref row, state, traceDict);
699            cur = AggregateApply(cur,
700              s => double.NaN,
701              v => Statistics.PopulationKurtosis(v));
702            TraceEvaluation(currentInstr, cur);
703            return cur;
704          }
705        case OpCodes.EuclideanDistance: {
706            var x1 = Evaluate(dataset, ref row, state, traceDict);
707            var x2 = Evaluate(dataset, ref row, state, traceDict);
708            var cur = AggregateMultipleApply(x1, x2,
709              (lhs, rhs) => ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 0.0),
710              (s1, s2) => s1 - s2,
711              (s1, v2) => Math.Sqrt((s1 - v2).PointwisePower(2).Sum()),
712              (v1, s2) => Math.Sqrt((v1 - s2).PointwisePower(2).Sum()),
713              (v1, v2) => Math.Sqrt((v1 - v2).PointwisePower(2).Sum()));
714            TraceEvaluation(currentInstr, cur);
715            return cur;
716          }
717        case OpCodes.Covariance: {
718            var x1 = Evaluate(dataset, ref row, state, traceDict);
719            var x2 = Evaluate(dataset, ref row, state, traceDict);
720            var cur = AggregateMultipleApply(x1, x2,
721              (lhs, rhs) => ApplyVectorLengthStrategy(DifferentVectorLengthStrategy, lhs, rhs, 0.0),
722              (s1, s2) => 0,
723              (s1, v2) => 0,
724              (v1, s2) => 0,
725              (v1, v2) => Statistics.PopulationCovariance(v1, v2));
726            TraceEvaluation(currentInstr, cur);
727            return cur;
728          }
729        case OpCodes.SubVector: {
730            var cur = Evaluate(dataset, ref row, state, traceDict);
731            return WindowedFunctionApply(cur, (WindowedSymbolTreeNode)currentInstr.dynamicNode,
732              s => s,
733              v => v);
734          }
735        case OpCodes.Variable: {
736            if (row < 0 || row >= dataset.Rows) return EvaluationResult.NaN;
737            var variableTreeNode = (VariableTreeNode)currentInstr.dynamicNode;
738            if (currentInstr.data is IList<double> doubleList) {
739              var cur = new EvaluationResult(doubleList[row] * variableTreeNode.Weight);
740              TraceEvaluation(currentInstr, cur);
741              return cur;
742            }
743            if (currentInstr.data is IList<DoubleVector> doubleVectorList) {
744              var cur = new EvaluationResult(doubleVectorList[row] * variableTreeNode.Weight);
745              TraceEvaluation(currentInstr, cur);
746              return cur;
747            }
748            throw new NotSupportedException($"Unsupported type of variable: {currentInstr.data.GetType().GetPrettyName()}");
749          }
750        case OpCodes.BinaryFactorVariable: {
751            if (row < 0 || row >= dataset.Rows) return EvaluationResult.NaN;
752            var factorVarTreeNode = currentInstr.dynamicNode as BinaryFactorVariableTreeNode;
753            var cur = new EvaluationResult(((IList<string>)currentInstr.data)[row] == factorVarTreeNode.VariableValue ? factorVarTreeNode.Weight : 0);
754            TraceEvaluation(currentInstr, cur);
755            return cur;
756          }
757        case OpCodes.FactorVariable: {
758            if (row < 0 || row >= dataset.Rows) return EvaluationResult.NaN;
759            var factorVarTreeNode = currentInstr.dynamicNode as FactorVariableTreeNode;
760            var cur = new EvaluationResult(factorVarTreeNode.GetValue(((IList<string>)currentInstr.data)[row]));
761            TraceEvaluation(currentInstr, cur);
762            return cur;
763          }
764        case OpCodes.Constant: {
765            var constTreeNode = (ConstantTreeNode)currentInstr.dynamicNode;
766            var cur = new EvaluationResult(constTreeNode.Value);
767            TraceEvaluation(currentInstr, cur);
768            return cur;
769          }
770
771        #region Time Series Symbols
772        case OpCodes.Median: {
773            var cur = Evaluate(dataset, ref row, state, traceDict);
774            cur = AggregateApply(cur,
775              s => s,
776              v => Statistics.Median(v));
777            TraceEvaluation(currentInstr, cur);
778            return cur;
779          }
780        case OpCodes.Quantile: {
781            var cur = Evaluate(dataset, ref row, state, traceDict);
782            var q = Evaluate(dataset, ref row, state, traceDict);
783            cur = AggregateApply(cur,
784              s => s,
785              v => Statistics.Quantile(v, q.Scalar));
786            TraceEvaluation(currentInstr, cur);
787            return cur;
788          }
789
790        case OpCodes.AbsoluteEnergy: {
791            var cur = Evaluate(dataset, ref row, state, traceDict);
792            cur = AggregateApply(cur,
793              s => s * s,
794              v => v.PointwisePower(2.0).Sum());
795            TraceEvaluation(currentInstr, cur);
796            return cur;
797          }
798
799        case OpCodes.BinnedEntropy: {
800            var cur = Evaluate(dataset, ref row, state, traceDict);
801            var m = Evaluate(dataset, ref row, state, traceDict);
802            cur = AggregateApply(cur,
803              s => 0,
804              v => {
805                int bins = Math.Max((int)Math.Round(m.Scalar), 1);
806                double minValue = v.Minimum();
807                double maxValue = v.Maximum();
808                double intervalWidth = (maxValue - minValue) / bins;
809                int totalValues = v.Count;
810                double sum = 0;
811                for (int i = 0; i < Math.Max(bins, v.Count); i++) {
812                  double binMin = minValue * i;
813                  double binMax = binMin + intervalWidth;
814                  double countBin = v.Map(e => (e > binMin && e < binMax) ? 1.0 : 0.0).Sum();
815                  double percBin = countBin / totalValues;
816                  sum += percBin * Math.Log(percBin);
817                }
818
819                return sum;
820              });
821            TraceEvaluation(currentInstr, cur);
822            return cur;
823          }
824        case OpCodes.HasLargeStandardDeviation: {
825            var cur = Evaluate(dataset, ref row, state, traceDict);
826            cur = AggregateApply(cur,
827              s => 0,
828              v => Statistics.PopulationStandardDeviation(v) > (Statistics.Maximum(v) - Statistics.Minimum(v)) / 2 ? 1.0 : 0.0);
829            TraceEvaluation(currentInstr, cur);
830            return cur;
831          }
832        case OpCodes.HasVarianceLargerThanStd: {
833            var cur = Evaluate(dataset, ref row, state, traceDict);
834            cur = AggregateApply(cur,
835              s => 0,
836              v => Statistics.PopulationVariance(v) > Statistics.StandardDeviation(v) ? 1.0 : 0.0);
837            TraceEvaluation(currentInstr, cur);
838            return cur;
839          }
840        case OpCodes.IsSymmetricLooking: {
841            var cur = Evaluate(dataset, ref row, state, traceDict);
842            cur = AggregateApply(cur,
843              s => 0,
844              v => Math.Abs(Statistics.Mean(v) - Statistics.Median(v)) < (Statistics.Maximum(v) - Statistics.Minimum(v)) / 2 ? 1.0 : 0.0);
845            TraceEvaluation(currentInstr, cur);
846            return cur;
847          }
848        case OpCodes.NumberDataPointsAboveMean: {
849            var cur = Evaluate(dataset, ref row, state, traceDict);
850            cur = AggregateApply(cur,
851              s => 0,
852              v => {
853                double mean = Statistics.Mean(v);
854                return v.Map(e => e > mean ? 1.0 : 0.0).Sum();
855              });
856            TraceEvaluation(currentInstr, cur);
857            return cur;
858          }
859        case OpCodes.NumberDataPointsAboveMedian: {
860            var cur = Evaluate(dataset, ref row, state, traceDict);
861            cur = AggregateApply(cur,
862              s => 0,
863              v => {
864                double median = Statistics.Median(v);
865                return v.Map(e => e > median ? 1.0 : 0.0).Sum();
866              });
867            TraceEvaluation(currentInstr, cur);
868            return cur;
869          }
870        case OpCodes.NumberDataPointsBelowMean: {
871            var cur = Evaluate(dataset, ref row, state, traceDict);
872            cur = AggregateApply(cur,
873              s => 0,
874              v => {
875                double mean = Statistics.Mean(v);
876                return v.Map(e => e < mean ? 1.0 : 0.0).Sum();
877              });
878            TraceEvaluation(currentInstr, cur);
879            return cur;
880          }
881        case OpCodes.NumberDataPointsBelowMedian: {
882            var cur = Evaluate(dataset, ref row, state, traceDict);
883            cur = AggregateApply(cur,
884              s => 0,
885              v => {
886                double median = Statistics.Median(v);
887                return v.Map(e => e < median ? 1.0 : 0.0).Sum();
888              });
889            TraceEvaluation(currentInstr, cur);
890            return cur;
891          }
892
893        case OpCodes.ArimaModelCoefficients: {
894            var cur = Evaluate(dataset, ref row, state, traceDict);
895            var i = Evaluate(dataset, ref row, state, traceDict);
896            var k = Evaluate(dataset, ref row, state, traceDict);
897            cur = AggregateApply(cur,
898              s => 0,
899              v => throw new NotImplementedException(""));
900            TraceEvaluation(currentInstr, cur);
901            return cur;
902          }
903        case OpCodes.ContinuousWaveletTransformationCoefficients: {
904            var cur = Evaluate(dataset, ref row, state, traceDict);
905            var a = Evaluate(dataset, ref row, state, traceDict);
906            var b = Evaluate(dataset, ref row, state, traceDict);
907            cur = AggregateApply(cur,
908              s => 0,
909              v => throw new NotImplementedException(""));
910            TraceEvaluation(currentInstr, cur);
911            return cur;
912          }
913        case OpCodes.FastFourierTransformationCoefficient: {
914            var cur = Evaluate(dataset, ref row, state, traceDict);
915            var k = Evaluate(dataset, ref row, state, traceDict);
916            cur = AggregateApply(cur,
917              s => 0,
918              v => throw new NotImplementedException(""));
919            TraceEvaluation(currentInstr, cur);
920            return cur;
921          }
922        case OpCodes.FirstIndexMax: {
923            var cur = Evaluate(dataset, ref row, state, traceDict);
924            cur = AggregateApply(cur,
925              s => 0,
926              v => (double)v.MaximumIndex() / v.Count);
927            TraceEvaluation(currentInstr, cur);
928            return cur;
929          }
930        case OpCodes.FirstIndexMin: {
931            var cur = Evaluate(dataset, ref row, state, traceDict);
932            cur = AggregateApply(cur,
933              s => 0,
934              v => (double)v.MinimumIndex() / v.Count);
935            TraceEvaluation(currentInstr, cur);
936            return cur;
937          }
938        case OpCodes.LastIndexMax: {
939            var cur = Evaluate(dataset, ref row, state, traceDict);
940            cur = AggregateApply(cur,
941              s => 0,
942              v => (double)(v.Count - DoubleVector.Build.DenseOfEnumerable(v.Reverse()).MaximumIndex()) / v.Count);
943
944            TraceEvaluation(currentInstr, cur);
945            return cur;
946          }
947        case OpCodes.LastIndexMin: {
948            var cur = Evaluate(dataset, ref row, state, traceDict);
949            cur = AggregateApply(cur,
950              s => 0,
951              v => (double)(v.Count - DoubleVector.Build.DenseOfEnumerable(v.Reverse()).MinimumIndex()) / v.Count);
952            TraceEvaluation(currentInstr, cur);
953            return cur;
954          }
955        case OpCodes.LongestStrikeAboveMean: {
956            var cur = Evaluate(dataset, ref row, state, traceDict);
957            cur = AggregateApply(cur,
958              s => 0,
959              v => LongestStrikeAbove(v, Statistics.Mean(v)));
960            TraceEvaluation(currentInstr, cur);
961            return cur;
962          }
963        case OpCodes.LongestStrikeAboveMedian: {
964            var cur = Evaluate(dataset, ref row, state, traceDict);
965            cur = AggregateApply(cur,
966              s => 0,
967              v => LongestStrikeAbove(v, Statistics.Median(v)));
968            TraceEvaluation(currentInstr, cur);
969            return cur;
970          }
971        case OpCodes.LongestStrikeBelowMean: {
972            var cur = Evaluate(dataset, ref row, state, traceDict);
973            cur = AggregateApply(cur,
974              s => 0,
975              v => LongestStrikeBelow(v, Statistics.Mean(v)));
976            TraceEvaluation(currentInstr, cur);
977            return cur;
978          }
979        case OpCodes.LongestStrikeBelowMedian: {
980            var cur = Evaluate(dataset, ref row, state, traceDict);
981            cur = AggregateApply(cur,
982              s => 0,
983              v => LongestStrikeBelow(v, Statistics.Median(v)));
984            TraceEvaluation(currentInstr, cur);
985            return cur;
986          }
987        case OpCodes.LongestStrikePositive: {
988            var cur = Evaluate(dataset, ref row, state, traceDict);
989            cur = AggregateApply(cur,
990              s => 0,
991              v => LongestStrikeAbove(v, 0));
992            TraceEvaluation(currentInstr, cur);
993            return cur;
994          }
995        case OpCodes.LongestStrikeNegative: {
996            var cur = Evaluate(dataset, ref row, state, traceDict);
997            cur = AggregateApply(cur,
998              s => 0,
999              v => LongestStrikeAbove(v, 0));
1000            TraceEvaluation(currentInstr, cur);
1001            return cur;
1002          }
1003        case OpCodes.LongestStrikeZero: {
1004            var cur = Evaluate(dataset, ref row, state, traceDict);
1005            cur = AggregateApply(cur,
1006              s => 0,
1007              v => LongestStrikeEqual(v, 0));
1008            TraceEvaluation(currentInstr, cur);
1009            return cur;
1010          }
1011        case OpCodes.MeanAbsoluteChange: {
1012            var cur = Evaluate(dataset, ref row, state, traceDict);
1013            cur = AggregateApply(cur,
1014              s => 0,
1015              v => {
1016                double sum = 0.0;
1017                for (int i = 0; i < v.Count - 1; i++) {
1018                  sum += Math.Abs(v[i + 1] - v[i]);
1019                }
1020
1021                return sum / v.Count;
1022              });
1023            TraceEvaluation(currentInstr, cur);
1024            return cur;
1025          }
1026        case OpCodes.MeanAbsoluteChangeQuantiles: {
1027            var cur = Evaluate(dataset, ref row, state, traceDict);
1028            var ql = Evaluate(dataset, ref row, state, traceDict);
1029            var qu = Evaluate(dataset, ref row, state, traceDict);
1030            cur = AggregateApply(cur,
1031              s => 0,
1032              v => {
1033                var lowerBound = Statistics.Quantile(v, ql.Scalar);
1034                var upperBound = Statistics.Quantile(v, qu.Scalar);
1035                var inBounds = v.Select(e => e > lowerBound && e < upperBound).ToList();
1036                double sum = 0.0;
1037                int count = 0;
1038                for (int i = 0; i < v.Count - 1; i++) {
1039                  if (inBounds[i] && inBounds[i + 1]) {
1040                    sum += Math.Abs(v[i + 1] - v[i]);
1041                    count++;
1042                  }
1043                }
1044
1045                return sum / count;
1046              });
1047            TraceEvaluation(currentInstr, cur);
1048            return cur;
1049          }
1050        case OpCodes.MeanAutocorrelation: {
1051            var cur = Evaluate(dataset, ref row, state, traceDict);
1052            cur = AggregateApply(cur,
1053              s => 0,
1054              v => {
1055                double sum = 0.0;
1056                double mean = Statistics.Mean(v);
1057                for (int l = 0; l < v.Count; l++) {
1058                  for (int i = 0; i < v.Count - l; i++) {
1059                    sum += (v[i] - mean) * (v[i + l] - mean);
1060                  }
1061                }
1062
1063                return sum / (v.Count - 1) / Statistics.PopulationVariance(v);
1064              });
1065            TraceEvaluation(currentInstr, cur);
1066            return cur;
1067          }
1068        case OpCodes.LaggedAutocorrelation: {
1069            var cur = Evaluate(dataset, ref row, state, traceDict);
1070            var lVal = Evaluate(dataset, ref row, state, traceDict);
1071            cur = AggregateApply(cur,
1072              s => 0,
1073              v => {
1074                double sum = 0.0;
1075                int l = Math.Max((int)Math.Round(lVal.Scalar), 0);
1076                double mean = Statistics.Mean(v);
1077                for (int i = 0; i < v.Count - l; i++) {
1078                  sum += (v[i] - mean) * (v[i + l] - mean);
1079                }
1080
1081                return sum / Statistics.PopulationVariance(v);
1082              });
1083            TraceEvaluation(currentInstr, cur);
1084            return cur;
1085          }
1086        case OpCodes.MeanSecondDerivateCentral: {
1087            var cur = Evaluate(dataset, ref row, state, traceDict);
1088            cur = AggregateApply(cur,
1089              s => 0,
1090              v => {
1091                double sum = 0.0;
1092                for (int i = 1; i < v.Count - 1; i++) {
1093                  sum += (v[i - 1] - 2 * v[i] + v[i + 1]) / 2;
1094                }
1095
1096                return sum / (v.Count - 2);
1097              });
1098            TraceEvaluation(currentInstr, cur);
1099            return cur;
1100          }
1101        case OpCodes.NumberPeaksOfSize: {
1102            var cur = Evaluate(dataset, ref row, state, traceDict);
1103            var l = Evaluate(dataset, ref row, state, traceDict);
1104            cur = AggregateApply(cur,
1105              s => 0,
1106              v => CountNumberOfPeaks(v, l.Scalar));
1107            TraceEvaluation(currentInstr, cur);
1108            return cur;
1109          }
1110        case OpCodes.LargeNumberOfPeaks: {
1111            var cur = Evaluate(dataset, ref row, state, traceDict);
1112            var l = Evaluate(dataset, ref row, state, traceDict);
1113            var m = Evaluate(dataset, ref row, state, traceDict);
1114            cur = AggregateApply(cur,
1115              s => 0,
1116              v => CountNumberOfPeaks(v, l.Scalar) > m.Scalar ? 1.0 : 0.0);
1117            TraceEvaluation(currentInstr, cur);
1118            return cur;
1119          }
1120        case OpCodes.TimeReversalAsymmetryStatistic: {
1121            var cur = Evaluate(dataset, ref row, state, traceDict);
1122            var l = Evaluate(dataset, ref row, state, traceDict);
1123            cur = AggregateApply(cur,
1124              s => 0,
1125              v => {
1126                int lag = Math.Max((int)Math.Round(l.Scalar), 0);
1127                double sum = 0.0;
1128                for (int i = 0; i < v.Count - 2 * lag; i++) {
1129                  sum += Math.Pow(v[i + 2 * lag], 2) * v[i + lag] - v[i + lag] * Math.Pow(v[i], 2);
1130                }
1131
1132                return sum / (v.Count - 2 * lag);
1133              });
1134            TraceEvaluation(currentInstr, cur);
1135            return cur;
1136          }
1137        #endregion
1138
1139        default:
1140          throw new NotSupportedException($"Unsupported OpCode: {currentInstr.opCode}");
1141      }
1142    }
1143
1144    private static int LongestStrikeAbove(DoubleVector v, double threshold) {
1145      int longestStrike = 0, currentStrike = 0;
1146      for (int i = 0; i < v.Count; i++) {
1147        if (v[i] > threshold) {
1148          currentStrike++;
1149          longestStrike = Math.Max(longestStrike, currentStrike);
1150        } else
1151          currentStrike = 0;
1152      }
1153      return longestStrike;
1154    }
1155    private static int LongestStrikeBelow(DoubleVector v, double threshold) {
1156      int longestStrike = 0, currentStrike = 0;
1157      for (int i = 0; i < v.Count; i++) {
1158        if (v[i] < threshold) {
1159          currentStrike++;
1160          longestStrike = Math.Max(longestStrike, currentStrike);
1161        } else
1162          currentStrike = 0;
1163      }
1164      return longestStrike;
1165    }
1166
1167    private static int LongestStrikeEqual(DoubleVector v, double value, double epsilon = double.Epsilon) {
1168      int longestStrike = 0, currentStrike = 0;
1169      for (int i = 0; i < v.Count; i++) {
1170        if (v[i].IsAlmost(epsilon)) {
1171          currentStrike++;
1172          longestStrike = Math.Max(longestStrike, currentStrike);
1173        } else
1174          currentStrike = 0;
1175      }
1176      return longestStrike;
1177    }
1178    private static int CountNumberOfPeaks(DoubleVector v, double heightDifference) {
1179      int count = 0;
1180      for (int i = 0; i < v.Count; i++) {
1181        bool largerThanPrev = i == 0 || v[i] > v[i - 1] + heightDifference;
1182        bool largerThanNext = i == v.Count - 1 || v[i] > v[i + 1] + heightDifference;
1183        if (largerThanPrev && largerThanNext)
1184          count++;
1185      }
1186      return count;
1187    }
1188  }
1189}
Note: See TracBrowser for help on using the repository browser.