diff --git a/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/DynamicTimeWarping.kt b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/DynamicTimeWarping.kt index 9ab96998e..72737ae57 100644 --- a/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/DynamicTimeWarping.kt +++ b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/DynamicTimeWarping.kt @@ -14,17 +14,9 @@ import kotlin.math.abs /** - * Offset constants which will be used later. Added them for avoiding "magical numbers" problem. - */ -internal enum class DtwOffset { - LEFT, - BOTTOM, - DIAGONAL -} - -/** - * Public class to store result of method. Class contains total penalty cost for series alignment. - * Also, this class contains align matrix (which point of the first series matches to point of the other series). + * Stores a result of [dynamicTimeWarping]. The class contains: + * 1. [Total penalty cost][totalCost] for series alignment. + * 2. [Align matrix][alignMatrix] that describes which point of the first series matches to point of the other series. */ public data class DynamicTimeWarpingData( val totalCost : Double = 0.0, @@ -32,41 +24,21 @@ public data class DynamicTimeWarpingData( ) /** - * PathIndices class for better code perceptibility. - * Special fun moveOption represent offset for indices. Arguments of this function - * is flags for bottom, diagonal or left offsets respectively. - */ -internal data class PathIndices (var id_x: Int, var id_y: Int) { - fun moveOption (direction: DtwOffset) { - when(direction) { - DtwOffset.BOTTOM -> id_x-- - DtwOffset.DIAGONAL -> { - id_x-- - id_y-- - } - DtwOffset.LEFT -> id_y-- - } - } -} - -/** - * Final DTW method realization. Returns alignment matrix - * for two series comparing and penalty for this alignment. + * DTW method implementation. Returns alignment matrix for two series comparing and penalty for this alignment. */ @OptIn(PerformancePitfall::class) public fun DoubleFieldOpsND.dynamicTimeWarping(series1 : DoubleBuffer, series2 : DoubleBuffer) : DynamicTimeWarpingData { - var cost = 0.0 - var pathLength = 0 - // Special matrix of costs alignment for two series. - val costMatrix = structureND(ShapeND(series1.size, series2.size)) { - (row, col) -> abs(series1[row] - series2[col]) + // Create a special matrix of costs alignment for the two series. + val costMatrix = structureND(ShapeND(series1.size, series2.size)) { (row, col) -> + abs(series1[row] - series2[col]) } - // Formula: costMatrix[i, j] = euqlideanNorm(series1(i), series2(j)) + - // min(costMatrix[i - 1, j], costMatrix[i, j - 1], costMatrix[i - 1, j - 1]). - for ( (row, col) in costMatrix.indices) { + + // Initialise the cost matrix by formulas + // costMatrix[i, j] = euclideanNorm(series1(i), series2(j)) + + // min(costMatrix[i - 1, j], costMatrix[i, j - 1], costMatrix[i - 1, j - 1]). + for ((row, col) in costMatrix.indices) { costMatrix[row, col] += when { - // There is special cases for i = 0 or j = 0. - row == 0 && col == 0 -> 0.0 + row == 0 && col == 0 -> continue row == 0 -> costMatrix[row, col - 1] col == 0 -> costMatrix[row - 1, col] else -> minOf( @@ -76,33 +48,37 @@ public fun DoubleFieldOpsND.dynamicTimeWarping(series1 : DoubleBuffer, series2 : ) } } + // alignMatrix contains non-zero values at position where two points from series matches // Values are penalty for concatenation of current points. - val alignMatrix = structureND(ShapeND(series1.size, series2.size)) {(_, _) -> 0.0} - val indexes = PathIndices(series1.size - 1, series2.size - 1) + val alignMatrix = structureND(ShapeND(series1.size, series2.size)) { _ -> 0.0} + var index1 = series1.size - 1 + var index2 = series2.size - 1 + var cost = 0.0 + var pathLength = 0 - with(indexes) { - alignMatrix[id_x, id_y] = costMatrix[id_x, id_y] - cost += costMatrix[id_x, id_y] - pathLength++ - while (id_x != 0 || id_y != 0) { - when { - id_x == 0 || costMatrix[id_x, id_y] == costMatrix[id_x, id_y - 1] + abs(series1[id_x] - series2[id_y]) -> { - moveOption(DtwOffset.LEFT) - } - id_y == 0 || costMatrix[id_x, id_y] == costMatrix[id_x - 1, id_y] + abs(series1[id_x] - series2[id_y]) -> { - moveOption(DtwOffset.BOTTOM) - } - costMatrix[id_x, id_y] == costMatrix[id_x - 1, id_y - 1] + abs(series1[id_x] - series2[id_y]) -> { - moveOption(DtwOffset.DIAGONAL) - } + alignMatrix[index1, index2] = costMatrix[index1, index2] + cost += costMatrix[index1, index2] + pathLength++ + while (index1 != 0 || index2 != 0) { + when { + index1 == 0 || costMatrix[index1, index2] == costMatrix[index1, index2 - 1] + abs(series1[index1] - series2[index2]) -> { + index2-- + } + index2 == 0 || costMatrix[index1, index2] == costMatrix[index1 - 1, index2] + abs(series1[index1] - series2[index2]) -> { + index1-- + } + costMatrix[index1, index2] == costMatrix[index1 - 1, index2 - 1] + abs(series1[index1] - series2[index2]) -> { + index1-- + index2-- } - alignMatrix[id_x, id_y] = costMatrix[id_x, id_y] - cost += costMatrix[id_x, id_y] - pathLength++ } - cost /= pathLength + alignMatrix[index1, index2] = costMatrix[index1, index2] + cost += costMatrix[index1, index2] + pathLength++ } + cost /= pathLength + return DynamicTimeWarpingData(cost, alignMatrix) } diff --git a/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/DTWTest.kt b/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/DTWTest.kt index 292fa5da2..6ff107426 100644 --- a/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/DTWTest.kt +++ b/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/DTWTest.kt @@ -15,25 +15,21 @@ import kotlin.test.Test class DTWTest { @Test - fun someData() : Unit { - with(Double.algebra.bufferAlgebra.seriesAlgebra()) { - val firstSequence: DoubleArray = doubleArrayOf(0.0, 2.0, 3.0, 1.0, 3.0, 0.1, 0.0, 1.0) - val secondSequence: DoubleArray = doubleArrayOf(1.0, 0.0, 3.0, 0.0, 0.0, 3.0, 2.0, 0.0, 2.0) + fun someData() { + val firstSequence: DoubleArray = doubleArrayOf(0.0, 2.0, 3.0, 1.0, 3.0, 0.1, 0.0, 1.0) + val secondSequence: DoubleArray = doubleArrayOf(1.0, 0.0, 3.0, 0.0, 0.0, 3.0, 2.0, 0.0, 2.0) - val seriesOne = firstSequence.asBuffer() - val seriesTwo = secondSequence.asBuffer() + val seriesOne = firstSequence.asBuffer() + val seriesTwo = secondSequence.asBuffer() - val result = DoubleFieldOpsND.dynamicTimeWarping(seriesOne, seriesTwo) - println("Total penalty coefficient: ${result.totalCost}") - print("Alignment: ") - println(result.alignMatrix) - for ((i , j) in result.alignMatrix.indices) { - if (result.alignMatrix[i, j] > 0.0) { - print("[$i, $j] ") - } + val result = DoubleFieldOpsND.dynamicTimeWarping(seriesOne, seriesTwo) + println("Total penalty coefficient: ${result.totalCost}") + print("Alignment: ") + println(result.alignMatrix) + for ((i , j) in result.alignMatrix.indices) { + if (result.alignMatrix[i, j] > 0.0) { + print("[$i, $j] ") } } } } - -