diff --git a/scripts/builtin/hmm.dml b/scripts/builtin/hmm.dml new file mode 100644 index 00000000000..02fafa6c054 --- /dev/null +++ b/scripts/builtin/hmm.dml @@ -0,0 +1,258 @@ +#------------------------------------------------------------- +# +# Licensed to the Apache Software Foundation (ASF) under one +# or more contributor license agreements. See the NOTICE file +# distributed with this work for additional information +# regarding copyright ownership. The ASF licenses this file +# to you under the Apache License, Version 2.0 (the +# "License"); you may not use this file except in compliance +# with the License. You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, +# software distributed under the License is distributed on an +# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +# KIND, either express or implied. See the License for the +# specific language governing permissions and limitations +# under the License. +# +#------------------------------------------------------------- + +# This script implements the hidden markov model method +# INPUT: +# -------------------------------------------------------------------------------------------- +# X Output of hmm of last n timestep +# OUTPUT: +# -------------------------------------------------------------------------------------------- +# ip Start probabilities +# A Transition probabilities +# B Emission probabilities +# -------------------------------------------------------------------------------------------- + +m_hmm = function(Matrix[Double] X) return (Matrix[Double] ip, Matrix[Double] A, Matrix[Double] B) +{ + + #X should have the size of 1 * ncols + #should be transposed for the unique function + unique_X = matrix(X, rows=ncol(X), cols=1) + nr_outputs = unique_vals(unique_X) + + /* + Since nr of states if unknown, fit a hmm model for every total number of states + from 1 to some max_states or 1 to log(n_timesteps), depending on whichever is greater. + + Also the break condition is: If the likelihood decreases with increase in number of total states, + break and take the paremeters of the last iteration as the optimal one. + */ + + max_states = 6 + T = ncol(X) + + if (max_states < log(T)){ + max_states = ceil(log(T)) + } + + search = TRUE + nr_states = 2 + seed = 42 + + while (search) { + [A_new, B_new, ip_new, curr_ll, seed] = baum_welch(X, nr_states, nr_outputs, 20, seed) + seed = seed+1 + if (nr_states == 2) { + prev_ll = -1 + } + if (curr_ll < prev_ll) { + search = FALSE + } + else { + A = A_new + B = B_new + ip = ip_new + } + prev_ll = curr_ll + nr_states = nr_states+1 + + if (nr_states > max_states) { + search = FALSE + } + } +} + +col_normalize = function(Matrix[Double] m) return (Matrix[Double] m_norm) { + m_norm = m / rowSums(m) +} + +unique_vals = function (Matrix[Double] X) return (Integer l) { + # group and count duplicates + I = aggregate (target = X, groups = X[,1], fn = "count") + # select groups + res = removeEmpty (target = seq (1, max (X[,1])), margin = "rows", select = (I != 0)) + l = length(res) +} + +baum_welch = function( Matrix[Double] X, Integer nr_states, Integer nr_outputs, Integer IT, Integer seed) return (Matrix[Double] A, Matrix[Double] B, Matrix[Double] ip, Double likelihood, Integer seed) +{ + #initialize state transition and emmission matrices uniformly + A = col_normalize(matrix(1/nr_states, rows=nr_states, cols=nr_states) + rand(rows=nr_states, cols=nr_states, seed=seed)) + seed = seed+1 + B = col_normalize(matrix(1/nr_outputs,rows=nr_states, cols=nr_outputs) + rand(rows=nr_states, cols=nr_outputs, seed=seed)) + seed = seed+1 + ip = matrix(1/nr_states, rows=nr_states, cols=1) + rand(rows=nr_states, cols=1, seed=seed) + seed = seed+1 + + likelihood = 0 + T = ncol(X) + + for (i in 1:IT) { + alpha = forward(X, A, B, ip) + beta = backward(X, A, B) + + gamma = calculate_gamma(alpha, beta) + eta = calculate_eta(X, alpha, beta, A, B) + + ip = gamma[, 1] + A = calculate_A(eta, gamma, nr_states) + B = calculate_B(gamma, X, nr_states, nr_outputs) + + #likelhood using forward method + likelihood = sum(gamma[, T] * alpha[, T]) + } +} + +calculate_A = function(Matrix[Double] eta, Matrix[Double] gamma, Integer nr_states) return (Matrix[Double] A) { + A = rand(rows=nr_states, cols=nr_states) + T = ncol(eta) + + parfor (i in 1:nr_states) { + parfor (j in 1:nr_states) { + eta_id = (i-1) * nr_states + j + num_ij = sum(eta[eta_id, 1:T-1]) + den_ij = sum(gamma[i,1:T-1]) + A[i,j] = num_ij/den_ij + } + } +} + +calculate_B = function(Matrix[Double] gamma, Matrix[Double] X, Integer nr_states, Integer nr_outputs) return (Matrix[Double] B) { + B = rand(rows=nr_states, cols=nr_outputs) + + parfor (i in 1:nr_states){ + parfor (j in 1:nr_outputs) { + B[i, j] = sum((X == j) * gamma[i, ]) / sum(gamma[i, ]) + } + } +} + +/* + To work with 3d matrices, the values in first two dimensions of eta matrix + must be rolled the first dimension. The following indexing function maps the + index of rolled matrix(eta) to "normal" 2d matrices like alpha, beta, etc +*/ +index = function(Integer rolled_i, Integer nr_states) return(Double i, Double j) { + i = ceil(rolled_i / nr_states) #index starts at 1 + j = rolled_i - (nr_states * (i-1)) #i-1) * +} + +calculate_eta = function (Matrix[Double] X, Matrix[Double] alpha, Matrix[Double] beta, Matrix[Double] A, Matrix[Double] B) return (Matrix[Double] eta) +{ + nr_states = nrow(alpha) + T = ncol(alpha) + + /* + eta a (nr_states * nr_states) * T matrix with cell (i, t) being probability of + the state being at i at timestep j given the observed output + */ + + tot_transitions = nr_states * nr_states + eta = rand(rows=tot_transitions, cols=T-1) + + /* + The transitions will be indiced as such:transition 1-N will represent + transition from state 1 to 1, 2 upto state N. transition N+1-2N will represent + transitions from 2 to 1, 2 upto state N + */ + + parfor (trans_id in 1:tot_transitions) { + [i, j] = index(trans_id, nr_states) + for (t in 1:(T-1)) { + #indices for alpha and beta + ot1 = output_t(X, t+1) + num_ij = alpha[i, t] * A[i, j] * beta[j, t+1] * B[j, ot1] + den = matrix(0, rows=nr_states, cols=1) + + parfor (k in 1:nr_states) { + den[k, 1] = sum(alpha[k, t] * t(A[k, ]) * beta[, t+1] * B[,ot1]) + } + s = num_ij / sum(den) + eta[trans_id, t] = as.scalar(s[1,1]) + } + } +} + +calculate_gamma = function (Matrix[Double] alpha, Matrix[Double] beta) return (Matrix[Double] gamma) +{ + /* + gamma a nr_state * T matrix with cell (i, t) being probability of + the state being at i at timestep j given the observed output + */ + nr_states = nrow(alpha) + T = ncol(alpha) + + gamma = rand(rows=nr_states, cols=T) + + parfor (i in 1:nrow(alpha)) { + parfor (t in 1:ncol(alpha)) { + num_it = alpha[i, t] * beta[i, t] + den_it = sum(alpha[,t] * beta[,t]) + gamma[i, t] = num_it / den_it + } + } +} + +output_t = function (Matrix[Double] X, Integer t) return (Integer o) { + o = as.integer(as.scalar(X[1, t])) +} + +forward = function (Matrix[Double] X, Matrix[Double] A, Matrix[Double] B, Matrix[Double] ip) return (Matrix[Double] alpha) +{ + /* + alpha a matrix of size nr_states * T with a cell i,t being probability of + the state being at state i at timestep j and the outputs till timestep j + */ + + T = ncol(X) + nr_states = nrow(A) + alpha = matrix(0, rows=nr_states, cols=T) + o_1 = output_t(X, 1) + alpha[ ,1] = ip * B[ ,o_1] + + for (t in 2:T) { + ot = output_t(X, t) + for (i in 1:nr_states) { + alpha[i, t] = B[i, ot] * sum(alpha[,t-1]* A[ ,i]) + } + } +} + +backward = function (Matrix[Double] X, Matrix[Double] A, Matrix[Double] B) return (Matrix[Double] beta) +{ + /* + alpha a matrix of size nr_states * T with a cell i,t being probability of + the model producing outputs (o_t+1,..., o_T) given that the model is + at state i at time t + */ + + T = ncol(X) + nr_states = nrow(A) + beta = matrix(1/T, rows=nr_states, cols=T) + beta[,T] = matrix(1/T, rows=nr_states, cols=1) + + for (t in (T-1):1) { + ot = output_t(X, t) + for (i in 1:nr_states) { + beta[i, t] = sum(beta[, t+1] * t(A[i, ]) * B[ , ot]) + } + } +} \ No newline at end of file diff --git a/scripts/builtin/hmmPredict.dml b/scripts/builtin/hmmPredict.dml new file mode 100644 index 00000000000..dd838d6d8fe --- /dev/null +++ b/scripts/builtin/hmmPredict.dml @@ -0,0 +1,76 @@ +#------------------------------------------------------------- +# +# Licensed to the Apache Software Foundation (ASF) under one +# or more contributor license agreements. See the NOTICE file +# distributed with this work for additional information +# regarding copyright ownership. The ASF licenses this file +# to you under the Apache License, Version 2.0 (the +# "License"); you may not use this file except in compliance +# with the License. You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, +# software distributed under the License is distributed on an +# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +# KIND, either express or implied. See the License for the +# specific language governing permissions and limitations +# under the License. +# +#------------------------------------------------------------- + +# This script implements the hidden markov model method +# INPUT: +# -------------------------------------------------------------------------------------------- +# X Output of hmm of last n timestep +# k Number of timesteps to predict the output for +# OUTPUT: +# -------------------------------------------------------------------------------------------- +# outputs Output of hmm for k timesteps +# -------------------------------------------------------------------------------------------- + +m_hmmPredict = function (Matrix[Double] X, Integer k) return (Matrix[Double] outputs) +{ + /* + Only discrete model implemented, need to discretize the cotinous outputs + to work with the continous model + */ + [ip, A, B] = hmm(X) + /* + Also only discrete model implemented, need to discretize the cotinous outputs + to work with the continous model + */ + seed = 42 + init_state = choose_value(ip, seed) + + states = matrix(0, rows=k, cols=1) + + states[1, 1] = init_state + s = init_state + for (i in 2:k) { + seed=seed +1 + next_state = choose_value(t(A[s, ]), seed) + states[i, 1] = next_state + s = next_state + } + + outputs = matrix(0, rows=nrow(states), cols=1) + for (i in 1:nrow(states)) { + seed = seed+1 + s = as.scalar(states[i, 1]) + outputs[i, 1] = choose_value(t(B[s, ]), seed) + } +} + +choose_value = function (Matrix[Double] prob_dist, Integer seed) return (Double val) +{ + # for selecting output or states given probability distribution + + # choosing using cumulative distribution + cum_prob = cumsum(prob_dist) + + R = rand(rows=1, cols=1, min=0, max=1, seed=seed) + + # adding 1 to colsum because indexing starts from 1 + val = as.scalar(colSums(cum_prob <= R) + 1) +} diff --git a/src/main/java/org/apache/sysds/common/Builtins.java b/src/main/java/org/apache/sysds/common/Builtins.java index efb260360e7..821c2f6bbb1 100644 --- a/src/main/java/org/apache/sysds/common/Builtins.java +++ b/src/main/java/org/apache/sysds/common/Builtins.java @@ -151,6 +151,8 @@ public enum Builtins { GNMF("gnmf", true), GRID_SEARCH("gridSearch", true), TOPK_CLEANING("topk_cleaning", true), + HMM("hmm", true), + HMMPREDICT("hmmPredict", true), HOSPITAL_RESIDENCY_MATCH("hospitalResidencyMatch", true), HYPERBAND("hyperband", true), IFELSE("ifelse", false), diff --git a/src/test/java/org/apache/sysds/test/functions/builtin/part1/BuiltinHMMPredictTest.java b/src/test/java/org/apache/sysds/test/functions/builtin/part1/BuiltinHMMPredictTest.java new file mode 100644 index 00000000000..feb1174abd6 --- /dev/null +++ b/src/test/java/org/apache/sysds/test/functions/builtin/part1/BuiltinHMMPredictTest.java @@ -0,0 +1,101 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one + * or more contributor license agreements. See the NOTICE file + * distributed with this work for additional information + * regarding copyright ownership. The ASF licenses this file + * to you under the Apache License, Version 2.0 (the + * "License"); you may not use this file except in compliance + * with the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, + * software distributed under the License is distributed on an + * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY + * KIND, either express or implied. See the License for the + * specific language governing permissions and limitations + * under the License. + */ + +package org.apache.sysds.test.functions.builtin.part1; + +import java.util.HashMap; +import java.util.Random; +import org.apache.sysds.api.DMLScript; +import org.apache.sysds.common.Types; +import org.apache.sysds.common.Types.ExecType; +import org.apache.sysds.hops.OptimizerUtils; +import org.apache.sysds.runtime.matrix.data.MatrixValue; +import org.apache.sysds.test.AutomatedTestBase; +import org.apache.sysds.test.TestConfiguration; +import org.apache.sysds.test.TestUtils; +import org.junit.Test; + +@net.jcip.annotations.NotThreadSafe +public class BuiltinHMMPredictTest extends AutomatedTestBase +{ + private final static String TEST_NAME = "hmmPredict"; + private final static String TEST_DIR = "functions/builtin/"; + private static final String TEST_CLASS_DIR = TEST_DIR + BuiltinHMMPredictTest.class.getSimpleName() + "/"; + private final static double eps = 1e-3; + + @Override + public void setUp() { + TestUtils.clearAssertionInformation(); + addTestConfiguration(TEST_NAME,new TestConfiguration(TEST_CLASS_DIR, TEST_NAME,new String[]{"C"})); + } + + @Test + public void testHMMPredictCP() { + runHMMPredictTest(ExecType.CP); + } + + @Test + public void testHMMPredictSpark() { + runHMMPredictTest(ExecType.SPARK); + } + + private void runHMMPredictTest(ExecType instType) + { + Types.ExecMode platformOld = setExecMode(instType); + try + { + loadTestConfiguration(getTestConfiguration(TEST_NAME)); + String HOME = SCRIPT_DIR + TEST_DIR; + + fullDMLScriptName = HOME + TEST_NAME + ".dml"; + programArgs = new String[]{ + "-nvargs", "X=" + input("X"), "k=" + input("k"), "outputs=" + output("outputs")}; + + //Output generated from https://en.wikipedia.org/wiki/Hidden_Markov_model#Weather_guessing_game + double[][] X = {{2, 1, 2, 1, 2, 2, 3, 3, 1, 1, 1, 1, 2, 2, 3, 2, 1, 2, 3, 2, + 1, 1, 1, 3, 2, 2, 2, 2, 3, 3, 1, 2, 3, 1, 3, 3, 2, 3, 1, 2, 1, + 1, 2, 3, 1, 1, 1, 3, 2, 3, 2, 1, 1, 3, 3, 3, 1, 2, 3, 3, 1, 3, + 3, 3, 2, 2, 2, 3, 3, 3, 3, 3, 1, 1, 3, 2, 2, 3, 2, 2, 1, 3, 2, + 1, 3, 2, 2, 2, 1, 2, 1, 2, 2, 1, 3, 2, 2, 2, 3, 2}}; + int k = 10; + + writeInputMatrixWithMTD("X", X, true); + runTest(true, false, null, -1); + + HashMap test_output = readDMLMatrixFromOutputDir("outputs"); + HashMap expected_output = new HashMap<>(); + + expected_output.put(new MatrixValue.CellIndex(1, 1), 1.0); + expected_output.put(new MatrixValue.CellIndex(2, 1), 3.0); + expected_output.put(new MatrixValue.CellIndex(3, 1), 3.0); + expected_output.put(new MatrixValue.CellIndex(4, 1), 2.0); + expected_output.put(new MatrixValue.CellIndex(5, 1), 1.0); + expected_output.put(new MatrixValue.CellIndex(6, 1), 2.0); + expected_output.put(new MatrixValue.CellIndex(7, 1), 3.0); + expected_output.put(new MatrixValue.CellIndex(8, 1), 1.0); + expected_output.put(new MatrixValue.CellIndex(9, 1), 3.0); + expected_output.put(new MatrixValue.CellIndex(10, 1), 1.0); + + TestUtils.compareMatrices(expected_output, test_output, eps, "Expected-DML", "Actual-DML"); + } + finally { + rtplatform = platformOld; + } + } +} diff --git a/src/test/java/org/apache/sysds/test/functions/builtin/part1/BuiltinHMMTest.java b/src/test/java/org/apache/sysds/test/functions/builtin/part1/BuiltinHMMTest.java new file mode 100644 index 00000000000..31766c45031 --- /dev/null +++ b/src/test/java/org/apache/sysds/test/functions/builtin/part1/BuiltinHMMTest.java @@ -0,0 +1,121 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one + * or more contributor license agreements. See the NOTICE file + * distributed with this work for additional information + * regarding copyright ownership. The ASF licenses this file + * to you under the Apache License, Version 2.0 (the + * "License"); you may not use this file except in compliance + * with the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, + * software distributed under the License is distributed on an + * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY + * KIND, either express or implied. See the License for the + * specific language governing permissions and limitations + * under the License. + + */ + + +package org.apache.sysds.test.functions.builtin.part1; + + +import java.util.HashMap; +import java.util.Random; +import org.apache.sysds.api.DMLScript; +import org.apache.sysds.common.Types; +import org.apache.sysds.common.Types.ExecType; +import org.apache.sysds.hops.OptimizerUtils; +import org.apache.sysds.runtime.matrix.data.MatrixValue; +import org.apache.sysds.test.AutomatedTestBase; +import org.apache.sysds.test.TestConfiguration; +import org.apache.sysds.test.TestUtils; +import org.junit.Test; + +@net.jcip.annotations.NotThreadSafe +public class BuiltinHMMTest extends AutomatedTestBase +{ + private final static String TEST_NAME = "hmm"; + private final static String TEST_DIR = "functions/builtin/"; + private static final String TEST_CLASS_DIR = TEST_DIR + BuiltinHMMPredictTest.class.getSimpleName() + "/"; + private final static double eps = 1e-3; + + @Override + public void setUp() { + TestUtils.clearAssertionInformation(); + addTestConfiguration(TEST_NAME, new TestConfiguration(TEST_CLASS_DIR, TEST_NAME, new String[]{"C"})); + } + + @Test + public void testHMMPredictCP() { + runHMMPredictTest(ExecType.CP); + } + + @Test + public void testHMMPredictSpark() { + runHMMPredictTest(ExecType.SPARK); + } + + private void runHMMPredictTest(ExecType instType) + { + Types.ExecMode platformOld = setExecMode(instType); + try { + loadTestConfiguration(getTestConfiguration(TEST_NAME)); + + String HOME = SCRIPT_DIR + TEST_DIR; + fullDMLScriptName = HOME + TEST_NAME + ".dml"; + programArgs = new String[]{ + "-nvargs", "X=" + input("X"), "ip=" + output("ip"), "A=" + output("A"), "B=" + output("B")}; + + //Output generated from https://en.wikipedia.org/wiki/Hidden_Markov_model#Weather_guessing_game + double[][] X = {{2, 1, 2, 1, 2, 2, 3, 3, 1, 1, 1, 1, 2, 2, 3, 2, 1, 2, 3, 2, + 1, 1, 1, 3, 2, 2, 2, 2, 3, 3, 1, 2, 3, 1, 3, 3, 2, 3, 1, 2, 1, + 1, 2, 3, 1, 1, 1, 3, 2, 3, 2, 1, 1, 3, 3, 3, 1, 2, 3, 3, 1, 3, + 3, 3, 2, 2, 2, 3, 3, 3, 3, 3, 1, 1, 3, 2, 2, 3, 2, 2, 1, 3, 2, + 1, 3, 2, 2, 2, 1, 2, 1, 2, 2, 1, 3, 2, 2, 2, 3, 2}}; + + writeInputMatrixWithMTD("X", X, true); + runTest(true, false, null, -1); + + HashMap out_ip = readDMLMatrixFromOutputDir("ip"); + HashMap out_A = readDMLMatrixFromOutputDir("A"); + HashMap out_B = readDMLMatrixFromOutputDir("B"); + HashMap expected_ip = new HashMap<>(); + HashMap expected_A = new HashMap<>(); + HashMap expected_B = new HashMap<>(); + + expected_ip.put(new MatrixValue.CellIndex(1, 1), 0.731); + expected_ip.put(new MatrixValue.CellIndex(2, 1), 0.000); + expected_ip.put(new MatrixValue.CellIndex(3, 1), 0.269); + + expected_A.put(new MatrixValue.CellIndex(1, 1), 0.194); + expected_A.put(new MatrixValue.CellIndex(1, 2), 0.421); + expected_A.put(new MatrixValue.CellIndex(1, 3), 0.424); + expected_A.put(new MatrixValue.CellIndex(2, 1), 0.189); + expected_A.put(new MatrixValue.CellIndex(2, 2), 0.432); + expected_A.put(new MatrixValue.CellIndex(2, 3), 0.383); + expected_A.put(new MatrixValue.CellIndex(3, 1), 0.374); + expected_A.put(new MatrixValue.CellIndex(3, 2), 0.334); + expected_A.put(new MatrixValue.CellIndex(3, 3), 0.257); + + expected_B.put(new MatrixValue.CellIndex(1, 1), 0.383); + expected_B.put(new MatrixValue.CellIndex(1, 2), 0.609); + expected_B.put(new MatrixValue.CellIndex(1, 3), 0.007); + expected_B.put(new MatrixValue.CellIndex(2, 1), 0.152); + expected_B.put(new MatrixValue.CellIndex(2, 2), 0.012); + expected_B.put(new MatrixValue.CellIndex(2, 3), 0.837); + expected_B.put(new MatrixValue.CellIndex(3, 1), 0.373); + expected_B.put(new MatrixValue.CellIndex(3, 2), 0.615); + expected_B.put(new MatrixValue.CellIndex(3, 3), 0.012); + + TestUtils.compareMatrices(expected_ip, out_ip, eps, "Expected-DML", "Actual-DML"); + TestUtils.compareMatrices(expected_A, out_A, eps, "Expected-DML", "Actual-DML"); + TestUtils.compareMatrices(expected_B, out_B, eps, "Expected-DML", "Actual-DML"); + } + finally { + rtplatform = platformOld; + } + } +} diff --git a/src/test/scripts/functions/builtin/hmmPredictTest.dml b/src/test/scripts/functions/builtin/hmmPredictTest.dml new file mode 100644 index 00000000000..271be550317 --- /dev/null +++ b/src/test/scripts/functions/builtin/hmmPredictTest.dml @@ -0,0 +1,25 @@ +#------------------------------------------------------------- +# +# Licensed to the Apache Software Foundation (ASF) under one +# or more contributor license agreements. See the NOTICE file +# distributed with this work for additional information +# regarding copyright ownership. The ASF licenses this file +# to you under the Apache License, Version 2.0 (the +# "License"); you may not use this file except in compliance +# with the License. You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, +# software distributed under the License is distributed on an +# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +# KIND, either express or implied. See the License for the +# specific language governing permissions and limitations +# under the License. +# +#------------------------------------------------------------- + +X = read($1); +k = read($2); +Y = hmmPredict(X, k); +write(Y); diff --git a/src/test/scripts/functions/builtin/hmmTest.dml b/src/test/scripts/functions/builtin/hmmTest.dml new file mode 100644 index 00000000000..e69de29bb2d