diff --git a/HISTORY.txt b/HISTORY.txt
index a5155d7..610ffed 100644
--- a/HISTORY.txt
+++ b/HISTORY.txt
@@ -2,6 +2,11 @@ Changelog
 ==========
 
 
+15.0.0 (2026-08-14)
+-------------------
+
+* Initial release for DSS 15.0.0
+
 14.7.3 (2026-08-03)
 -------------------
 
diff --git a/dataikuscoring/algorithms/__init__.py b/dataikuscoring/algorithms/__init__.py
index 2872508..1572aec 100644
--- a/dataikuscoring/algorithms/__init__.py
+++ b/dataikuscoring/algorithms/__init__.py
@@ -3,6 +3,9 @@
 from .forest_regressor import ForestRegressor
 from .gradient_boosting_classifier import GradientBoostingClassifier
 from .gradient_boosting_regressor import GradientBoostingRegressor
+from .isolation_forest import IsolationForest
+from .kmeans import KMeans
+from .kmeans import MiniBatchKMeans
 from .linear_regression import LinearRegressor
 from .logistic import LogisticRegressionClassifier
 from .mlp_classifier import MLPClassifer
@@ -15,6 +18,9 @@
     "FOREST_REGRESSOR": ForestRegressor,
     "GRADIENT_BOOSTING_CLASSIFIER": GradientBoostingClassifier,
     "GRADIENT_BOOSTING_REGRESSOR": GradientBoostingRegressor,
+    "ISOLATION_FOREST": IsolationForest,
+    "KMEANS": KMeans,
+    "MINIBATCH_KMEANS": MiniBatchKMeans,
     "LINEAR": LinearRegressor,
     "LOGISTIC": LogisticRegressionClassifier,
     "MLP_REGRESSOR": MLPRegressor,
diff --git a/dataikuscoring/algorithms/common.py b/dataikuscoring/algorithms/common.py
index 25b9c8c..e9dfa9e 100644
--- a/dataikuscoring/algorithms/common.py
+++ b/dataikuscoring/algorithms/common.py
@@ -18,3 +18,18 @@ def __init__(self, model_parameters):
     def predict(self, X):
         """Predict target vector from a 2D numpy array input X"""
         raise NotImplementedError
+
+
+class Clusterer:
+
+    def __init__(self, model_parameters):
+        """The content of the dss_pipeline_model.gz file"""
+        raise NotImplementedError
+
+    def predict(self, X):
+        """Predict the cluster index for each row of a 2D numpy array input X.
+
+        Anomaly-detection clusterers (e.g. Isolation Forest) additionally expose decision_function(X)
+        returning the per-row anomaly score.
+        """
+        raise NotImplementedError
diff --git a/dataikuscoring/algorithms/decision_tree_model.py b/dataikuscoring/algorithms/decision_tree_model.py
index 84b4e34..0202d6a 100644
--- a/dataikuscoring/algorithms/decision_tree_model.py
+++ b/dataikuscoring/algorithms/decision_tree_model.py
@@ -59,7 +59,7 @@ class Node:
 
     def __init__(self, feature_idx=None, threshold=np.nan, left_child=None, right_child=None, label=None,
                  is_leaf=None, missing_goes_left=None, missing_value=np.nan, split_kind=SPLIT_KIND_THRESHOLD,
-                 category_set=None):
+                 category_set=None, n_node_samples=None):
         self.label = label
         self.feature_idx = feature_idx
         self.threshold = threshold
@@ -70,6 +70,8 @@ def __init__(self, feature_idx=None, threshold=np.nan, left_child=None, right_ch
         self.missing_value = missing_value
         self.split_kind = split_kind
         self.category_set = None if category_set is None else frozenset(float(v) for v in category_set)
+        # only populated for Isolation Forest leaves (used by its path-length anomaly score); None otherwise
+        self.n_node_samples = n_node_samples
 
     def is_missing(self, data):
         if np.isnan(self.missing_value):
@@ -136,9 +138,14 @@ def init_tree(self, model_parameters):
         missing_value = model_parameters.get("missing_value", np.nan)
 
         convert_threshold = np.float32 if self.variant == "XGBOOST" else np.float64
+        # n_node_samples is present only for Isolation Forest trees; aligned with leaf_id when present
+        leaf_n_node_samples = model_parameters.get("n_node_samples")
+        if leaf_n_node_samples is None or len(leaf_n_node_samples) == 0:
+            leaf_n_node_samples = [None] * len(model_parameters["leaf_id"])
         leaves = {
-            leaf_id: Node(label=label, is_leaf=True, missing_value=missing_value) for leaf_id, label in zip(
-                model_parameters["leaf_id"], model_parameters["label"])
+            leaf_id: Node(label=label, is_leaf=True, missing_value=missing_value, n_node_samples=n_node_samples)
+            for leaf_id, label, n_node_samples in zip(
+                model_parameters["leaf_id"], model_parameters["label"], leaf_n_node_samples)
         }
 
         missing = model_parameters.get("missing")
diff --git a/dataikuscoring/algorithms/isolation_forest.py b/dataikuscoring/algorithms/isolation_forest.py
new file mode 100644
index 0000000..247a5e0
--- /dev/null
+++ b/dataikuscoring/algorithms/isolation_forest.py
@@ -0,0 +1,73 @@
+import numpy as np
+
+from .common import Clusterer
+from .decision_tree_model import DecisionTreeModel, SPLIT_KIND_CATEGORY_SET
+
+
+def _average_path_length(n):
+    """c(n): expected path length of an unsuccessful search in a binary search tree (sklearn _average_path_length)."""
+    if n <= 1:
+        return 0.0
+    if n == 2:
+        return 1.0
+    return 2.0 * (np.log(n - 1) + np.euler_gamma) - 2.0 * (n - 1) / n
+
+
+def _path_depth(root, data):
+    """Walk a tree to its leaf using sklearn routing (<=, missing-goes-left), returning (edge count, leaf)."""
+    current = root
+    depth = 0
+    while not current.is_leaf:
+        if current.is_missing(data):
+            current = current.left_child if current.missing_goes_left else current.right_child
+        elif current.split_kind == SPLIT_KIND_CATEGORY_SET:
+            current = current.left_child if current.has_category(data) else current.right_child
+        elif data[current.feature_idx] <= current.threshold:
+            current = current.left_child
+        else:
+            current = current.right_child
+        depth += 1
+    return depth, current
+
+
+class IsolationForest(Clusterer):
+    """Anomaly-detection clusterer reproducing sklearn IsolationForest.
+
+    decision_function(X) returns the anomaly score (the doctor's anomaly_score): per tree the isolation path
+    is leaf_depth + c(n_node_samples@leaf), averaged over the trees, normalised by c(max_samples) into
+    2^(-mean/c(psi)); the returned score is -raw - offset (anomaly when < 0).
+
+    predict(X) returns the cluster index: 1 (anomaly) when the score is < 0, else 0 (regular) -- matching the
+    doctor's DkuIsolationForest.predict. The model layer maps these indices to the cluster names.
+    """
+
+    def __init__(self, model_parameters):
+        self.trees = [DecisionTreeModel(tree_params) for tree_params in model_parameters["trees"]]
+        # feature subset each tree was scored on (sklearn scores tree i on X[:, estimators_features[i]])
+        self.estimators_features = model_parameters["estimators_features"]
+        self.max_samples = model_parameters["max_samples"]
+        self.offset = model_parameters["offset"]
+        self._normalizer = _average_path_length(self.max_samples)
+        self.feature_converter = self.trees[0].feature_converter
+
+    def predict(self, X):
+        # cluster index: 1 -> anomaly (score < 0), 0 -> regular
+        return [int(score < 0) for score in self._scores(X)]
+
+    def decision_function(self, X):
+        return self._scores(X)
+
+    def _scores(self, X):
+        return [self._score(data) for data in self.feature_converter(X)]
+
+    def _score(self, data):
+        total_path = 0.0
+        for tree, features in zip(self.trees, self.estimators_features):
+            depth, leaf = _path_depth(tree.root, data[features])
+            total_path += depth + _average_path_length(leaf.n_node_samples)
+        mean_path = total_path / len(self.trees)
+        raw_score = 2.0 ** (-mean_path / self._normalizer)
+        return -raw_score - self.offset
+
+    def __repr__(self):
+        return "IsolationForest(n_trees={})".format(len(self.trees))
diff --git a/dataikuscoring/algorithms/kmeans.py b/dataikuscoring/algorithms/kmeans.py
new file mode 100644
index 0000000..1954054
--- /dev/null
+++ b/dataikuscoring/algorithms/kmeans.py
@@ -0,0 +1,38 @@
+import numpy as np
+
+from .common import Clusterer
+
+
+class KMeans(Clusterer):
+    """Clusterer reproducing sklearn KMeans / MiniBatchKMeans.
+
+    predict(X) assigns each row to the nearest cluster center (argmin squared-euclidean distance) in the
+    preprocessed feature space and returns the cluster index; the model layer maps that index to the cluster
+    name ("cluster_0", "cluster_1", ...). Both algorithms are fully defined by cluster_centers_, so one scorer
+    covers both. Unlike Isolation Forest there is no anomaly score, so no decision_function is exposed.
+    """
+
+    def __init__(self, model_parameters):
+        # k x n_features; row j is the center of cluster j
+        self.cluster_centers = np.asarray(model_parameters["cluster_centers"], dtype=np.float64)
+        # ||c||^2 per center, precomputed for the distance expansion below
+        self._center_sq_norms = (self.cluster_centers ** 2).sum(axis=1)
+
+    def predict(self, X):
+        data = np.asarray(X, dtype=np.float64)
+        # nearest center by squared-euclidean distance via the expansion ||x - c||^2 = ||x||^2 - 2 x.c + ||c||^2.
+        # ||x||^2 is constant per row so it is dropped (it does not change the argmin); this keeps the cost
+        # O(n*k) in memory and lets the Java engine reproduce the exact same formula for cross-engine parity.
+        distances = -2.0 * data.dot(self.cluster_centers.T) + self._center_sq_norms
+        return [int(i) for i in np.argmin(distances, axis=1)]
+
+    def __repr__(self):
+        return "KMeans(n_clusters={})".format(len(self.cluster_centers))
+
+
+class MiniBatchKMeans(KMeans):
+    """MiniBatchKMeans scores identically to KMeans (nearest center); the distinct class just preserves the
+    model type through a serialize/reload round-trip."""
+
+    def __repr__(self):
+        return "MiniBatchKMeans(n_clusters={})".format(len(self.cluster_centers))
diff --git a/dataikuscoring/algorithms/logistic.py b/dataikuscoring/algorithms/logistic.py
index 6114dbd..10f6a76 100644
--- a/dataikuscoring/algorithms/logistic.py
+++ b/dataikuscoring/algorithms/logistic.py
@@ -36,18 +36,20 @@ def multinomial_probabilities(dec):
 
 
 def modified_huber_probabilities(dec):
-    p = 0.5 * (1 + np.minimum(1, np.maximum(-1, dec)))
+    dec = np.asarray(dec, dtype=float)
+    p = 0.5 * (1 + np.clip(dec, -1, 1))
 
-    if len(dec[0]) == 2:
+    if dec.shape[1] == 2:
         p[:, 0] = 1 - p[:, 1]
 
-    norms = np.linalg.norm(dec, axis=1)
+    # scikit-learn normalizes the per-class values by their sum (Zadrozny & Elkan);
+    # rows whose values are all ~0 get a uniform distribution.
+    sums = p.sum(axis=1)
+    all_zero = sums < 1e-15
+    p[all_zero] = 1.0 / dec.shape[1]
+    sums[all_zero] = 1.0
 
-    # scikit-learn puts equal probas in this case
-    indexes = np.where(norms < 1e-15)
-    p[indexes] = np.ones(len(dec)) * (1 / len(dec))
-
-    return p / norms
+    return p / sums[:, None]
 
 
 POLICIES = {
diff --git a/dataikuscoring/load.py b/dataikuscoring/load.py
index 12ecdbf..89d08eb 100644
--- a/dataikuscoring/load.py
+++ b/dataikuscoring/load.py
@@ -157,7 +157,15 @@ def load_resources_from_resource_folder(resources_folder):
     user_meta_filename = os.path.join(resources_folder, "user_meta.json")
     if os.path.isfile(user_meta_filename):
         with open(user_meta_filename) as f:
-            resources["threshold"] = json.load(f).get("activeClassifierThreshold", 0.5)
+            user_meta = json.load(f)
+        resources["threshold"] = user_meta.get("activeClassifierThreshold", 0.5)
+        # User cluster renames (intrinsic cluster name -> user-chosen name). Mirrors the DSS python and java
+        # engines (reg_scoring_recipe / Build.remapClusterNames) so optimized scoring returns the same
+        # cluster_labels after a rename. Empty/absent for non-clustering models.
+        resources["cluster_name_map"] = {
+            cluster_id: cluster_data["name"]
+            for cluster_id, cluster_data in user_meta.get("clusterMetas", {}).items()
+        }
 
     return resources
 
@@ -247,11 +255,17 @@ def create_model(resources):
             algorithm_name = "MLP_REGRESSOR"
         else:
             algorithm_name = "MLP_CLASSIFIER"
+    # Apply user cluster renames to the intrinsic cluster names, mirroring the DSS python/java engines.
+    # No-op for non-clustering models (cluster_name_map is empty).
+    classes = resources["meta"].get("classes")
+    cluster_name_map = resources.get("cluster_name_map")
+    if classes is not None and cluster_name_map:
+        classes = [cluster_name_map.get(name, name) for name in classes]
     parameters = {
         "prepare_input": PrepareInput(resources),
         "algorithm": ALGORITHMS[algorithm_name](dict({"missing_value": resources["missing_value"]}, **resources["model_parameters"])),
         "preprocessings": Preprocessings(resources),
-        "classes": resources["meta"].get("classes"),
+        "classes": classes,
         "calibration": Calibrator(resources),
         "drop_rows": DropRows(resources)
     }
diff --git a/dataikuscoring/mlflow/classification.py b/dataikuscoring/mlflow/classification.py
index 1c17f01..19704af 100644
--- a/dataikuscoring/mlflow/classification.py
+++ b/dataikuscoring/mlflow/classification.py
@@ -188,7 +188,7 @@ def mlflow_classification_predict_to_scoring_data(mlflow_model, imported_model_m
         logger.info("MLflow outputs integers, converting")
         preds = pd.Series(mlflow_raw_preds)
         pred_df = pd.DataFrame({"prediction": mlflow_raw_preds})
-        pred_df["prediction"].replace(int_to_label_map, inplace=True)
+        pred_df["prediction"] = pred_df["prediction"].astype(object).replace(int_to_label_map)
     elif (isinstance(first_value, float) or isinstance(first_value, np.floating)) and \
             imported_model_meta["predictionType"] == "BINARY_CLASSIFICATION":
         # only a column of floats ... probably prediction of class 1
@@ -223,7 +223,7 @@ def mlflow_classification_predict_to_scoring_data(mlflow_model, imported_model_m
             preds = (probas_one > threshold).astype(int)
             pred_df = pd.DataFrame({"prediction": preds})
             logger.debug("Computed pred df %s" % pred_df)
-            pred_df["prediction"].replace(int_to_label_map, inplace=True)
+            pred_df["prediction"] = pred_df["prediction"].astype(object).replace(int_to_label_map)
             logger.info("Computed cleanpred df %s" % pred_df["prediction"].dtype)
 
     try:
@@ -238,7 +238,7 @@ def mlflow_classification_predict_to_scoring_data(mlflow_model, imported_model_m
         exception_with_cause.__cause__ = e
         raise exception_with_cause
 
-    if probas is not None and np.isnan(probas.to_numpy()).any():
+    if probas is not None and np.isnan(probas.fillna(np.nan).to_numpy(dtype=float)).any():
         raise Exception("MLflow model predicted NaN probabilities")
 
     logger.debug("Final pred_df: %s " % pred_df)
diff --git a/dataikuscoring/mlflow/common.py b/dataikuscoring/mlflow/common.py
index 5401525..494f095 100644
--- a/dataikuscoring/mlflow/common.py
+++ b/dataikuscoring/mlflow/common.py
@@ -5,6 +5,7 @@
 import pandas as pd
 import numpy as np
 
+
 logger = logging.getLogger(__name__)
 
 class DisableMLflowTypeEnforcement(object):
diff --git a/dataikuscoring/models/__init__.py b/dataikuscoring/models/__init__.py
index 2628888..13e2004 100644
--- a/dataikuscoring/models/__init__.py
+++ b/dataikuscoring/models/__init__.py
@@ -1,13 +1,15 @@
 from .regression import RegressionModel
 from .binary import BinaryModel
 from .multiclass import MulticlassModel
+from .clustering import ClusteringModel
 from .partitioned import ClassificationPartitionedModel, RegressionPartitionedModel
 from .mlflow import MLflowModel
 
 MODELS = {
     "REGRESSION": RegressionModel,
     "BINARY_PROBABILISTIC": BinaryModel,
-    "MULTICLASS_PROBABILISTIC": MulticlassModel
+    "MULTICLASS_PROBABILISTIC": MulticlassModel,
+    "CLUSTERING": ClusteringModel
 }
 
 PARTITIONED_MODELS = {
diff --git a/dataikuscoring/models/clustering.py b/dataikuscoring/models/clustering.py
new file mode 100644
index 0000000..42113a3
--- /dev/null
+++ b/dataikuscoring/models/clustering.py
@@ -0,0 +1,44 @@
+import numpy as np
+
+from .common import PredictionModelMixin, check_input_data
+from .model import BaseModel
+
+
+class ClusteringModel(BaseModel, PredictionModelMixin):
+    """Clustering model for optimized scoring.
+
+    Like the in-DSS doctor (and dataikuscoring's classification models), the primary output of predict(X) is the
+    cluster label name (e.g. "regular"/"anomalies"), obtained by mapping the algorithm's cluster index through
+    the serialized cluster names (``classes``). For anomaly-detection clusterers the underlying anomaly score is
+    available separately via decision_function(X).
+    """
+
+    def __init__(self, prepare_input, preprocessings, algorithm, drop_rows, classes=None, **kwargs):
+        super(ClusteringModel, self).__init__(prepare_input, preprocessings, algorithm, drop_rows)
+        # cluster label names, e.g. ["regular", "anomalies"]; index i -> classes[i]. May be None for models
+        # serialized before cluster names were emitted, in which case predict() falls back to the raw index.
+        self.classes = classes
+
+    def _compute_predict(self, X):
+        X_processed, valid_rows_mask = self._compute_preprocessed(X)
+        y_pred = np.array([None] * len(X), dtype=object)
+        indices = self.algorithm.predict(X_processed)
+        if self.classes is not None:
+            y_pred[valid_rows_mask] = [self.classes[int(i)] for i in indices]
+        else:
+            y_pred[valid_rows_mask] = indices
+        return y_pred
+
+    def decision_function(self, X):
+        """Per-row anomaly score (available for anomaly-detection clusterers such as Isolation Forest)."""
+        if not hasattr(self.algorithm, "decision_function"):
+            raise NotImplementedError(
+                "decision_function is only available for anomaly-detection clustering models")
+        check_input_data(X)
+        X_processed, valid_rows_mask = self._compute_preprocessed(X)
+        scores = np.full(len(X), np.nan)
+        scores[valid_rows_mask] = self.algorithm.decision_function(X_processed)
+        return scores
+
+    def __repr__(self):
+        return "{} Clusterer".format(self.algorithm)
diff --git a/dataikuscoring/processors/prepare_input.py b/dataikuscoring/processors/prepare_input.py
index f881ee5..b4373aa 100644
--- a/dataikuscoring/processors/prepare_input.py
+++ b/dataikuscoring/processors/prepare_input.py
@@ -1,7 +1,7 @@
 import numpy as np
 
 from ..utils import IndexedMatrix
-
+from ..utils.pandas_compat import to_numpy
 
 class PrepareInput:
 
@@ -64,6 +64,10 @@ def process(self, X):
 
         get_column_copy = None
         if isinstance(X, (list, np.ndarray)):  # type is List[List] or numpy array or List[dict]
+            def is_missing(value):
+                is_nat = isinstance(value, (np.datetime64, np.timedelta64)) and np.isnat(value)
+                return value is None or value is np.nan or is_nat
+
             data = X[0]
             if not isinstance(data, dict):  # type is List[List] or numpy array
                 if isinstance(data, (list, np.ndarray)):
@@ -77,11 +81,20 @@ def process(self, X):
             else:  # Type is List[Dict] because we handled validation in check_input_data
                 get_column_copy =  lambda X, index_column, column: np.array([x.get(column) for x in X])
         else:  # Type is Dataframe because we handled validation in check_input_data
+            try:
+                import pandas as pd
+            except ImportError:
+                raise NotImplementedError("pandas is required when scoring a pandas.DataFrame")
+            is_missing = pd.isna
             missing_columns = set(self.mandatory_input_column_names).difference(set(X.columns))
             if len(missing_columns) > 0:
                 raise ValueError("Missing column(s) in input DataFrame: {}".format(
                     ",".join(missing_columns)))
-            get_column_copy = lambda X, index_column, column: X[column].values
+            def get_column_copy(X, index_column, column):
+                series = X[column]
+                if pd.api.types.is_extension_array_dtype(series) and pd.api.types.is_numeric_dtype(series):
+                    return to_numpy(series, dtype=np.float64, na_value=np.nan)
+                return to_numpy(series)
 
         # Fill the input columns data into the right matrices
         for (index_column, column) in enumerate(self.input_column_names):
@@ -95,9 +108,10 @@ def process(self, X):
                             X_non_numeric[:, column] = np.where(np.isnan(data), None, data.astype(str))
                         else :
                             X_non_numeric[:, column] = np.where(np.isnan(data), None, data)
-                    else:  # we have to convert empty string and nan to None
-                        X_non_numeric[:, column] = np.where(data.astype(str) == "", None, data)
-                        X_non_numeric[:, column] = np.where([x is np.nan for x in X_non_numeric[:, column]], None, X_non_numeric[:, column])
+                    else:  # Convert empty strings and missing values to None
+                        X_non_numeric[:, column] = [
+                            None if is_missing(value) or value == "" else value for value in data
+                        ]
                 else:
                     if np.issubdtype(data.dtype, np.number):  # if data is not numeric check empty string
                         X_numeric[:, column] = data
diff --git a/dataikuscoring/utils/pandas_compat.py b/dataikuscoring/utils/pandas_compat.py
new file mode 100644
index 0000000..e83a583
--- /dev/null
+++ b/dataikuscoring/utils/pandas_compat.py
@@ -0,0 +1,20 @@
+import numpy as np
+
+
+def to_numpy(obj, *args, **kwargs):
+    """
+    Convert a pandas object to a NumPy array across supported pandas versions.
+    Duplicated from dataiku.core.compat.pandas_compat.to_numpy to avoid dependency.
+    """
+    try:
+        return obj.to_numpy(*args, **kwargs)
+    except (AttributeError, TypeError):
+        na_value = kwargs.pop("na_value", None)
+        if na_value is not None:
+            obj = obj.fillna(na_value)
+
+        dtype = kwargs.pop("dtype", None)
+        copy = kwargs.pop("copy", False)
+        if args or dtype is not None or copy:
+            return np.array(obj.values, *args, dtype=dtype, copy=copy)
+        return obj.values
diff --git a/dataikuscoring/utils/prediction_result.py b/dataikuscoring/utils/prediction_result.py
index e011ed5..3ae5ef3 100644
--- a/dataikuscoring/utils/prediction_result.py
+++ b/dataikuscoring/utils/prediction_result.py
@@ -87,7 +87,13 @@ def concat(prediction_results):
 class PredictionResult(AbstractPredictionResult):
 
     def as_dataframe(self, for_json_serialization=False):
-        prediction_df = pd.DataFrame({PREDICTION: self.preds})
+        if self.multi_target_variables:
+            prediction_df = pd.DataFrame(
+                self.preds,
+                columns=["{}_{}".format(PREDICTION, target) for target in self.multi_target_variables]
+            )
+        else:
+            prediction_df = pd.DataFrame({PREDICTION: self.preds})
         if not self.has_prediction_intervals():
             return prediction_df
         intervals = self.prediction_intervals
@@ -99,13 +105,14 @@ def as_dataframe(self, for_json_serialization=False):
             prediction_df[PREDICTION_INTERVAL_UPPER] = intervals[:, 1]
         return prediction_df
 
-    def __init__(self, preds, prediction_intervals=None):
+    def __init__(self, preds, prediction_intervals=None, multi_target_variables=None):
         """
         :type preds: np.ndarray
         :type prediction_intervals: np.ndarray
         """
         self._preds = preds
         self._prediction_intervals = prediction_intervals
+        self.multi_target_variables = multi_target_variables
 
     @property
     def prediction_intervals(self):
@@ -146,7 +153,7 @@ def _concat(prediction_results):
             intervals_concat = np.concatenate([pr._prediction_intervals for pr in prediction_results])
         else:
             intervals_concat = None
-        return PredictionResult(preds_concat, intervals_concat)
+        return PredictionResult(preds_concat, intervals_concat, multi_target_variables=prediction_results[0].multi_target_variables)
 
 
 class ClassificationPredictionResult(AbstractPredictionResult):
diff --git a/setup.py b/setup.py
index 8e034b3..cde8d60 100644
--- a/setup.py
+++ b/setup.py
@@ -5,7 +5,7 @@
 long_description = (open('README.md').read() + '\n\n' +
                     open('HISTORY.txt').read())
 
-VERSION = "14.7.3"
+VERSION = "15.0.0"
 
 setuptools.setup(
     name='dataiku-scoring',