diff --git a/pyFTS/benchmarks/benchmarks.py b/pyFTS/benchmarks/benchmarks.py
index 5103f02..e43a9a6 100644
--- a/pyFTS/benchmarks/benchmarks.py
+++ b/pyFTS/benchmarks/benchmarks.py
@@ -19,6 +19,7 @@ from mpl_toolkits.mplot3d import Axes3D
 from pyFTS.probabilistic import ProbabilityDistribution
 from pyFTS.common import Transformations
 from pyFTS.models import song, chen, yu, ismailefendi, sadaei, hofts, pwfts, ifts, cheng, hwang
+from pyFTS.models.multivariate import mvfts, wmvfts, cmvfts
 from pyFTS.models.ensemble import ensemble
 from pyFTS.benchmarks import Measures, naive, arima, ResidualAnalysis, quantreg, knn
 from pyFTS.benchmarks import Util as bUtil
@@ -57,10 +58,16 @@ def get_benchmark_point_methods():
 def get_point_methods():
     """Return all FTS methods for point forecasting"""
     return [song.ConventionalFTS, chen.ConventionalFTS, yu.WeightedFTS, ismailefendi.ImprovedWeightedFTS,
-            cheng.TrendWeightedFTS, sadaei.ExponentialyWeightedFTS, hofts.HighOrderFTS, hwang.HighOrderFTS,
+            cheng.TrendWeightedFTS, sadaei.ExponentialyWeightedFTS,
+            hofts.HighOrderFTS, hofts.WeightedHighOrderFTS, hwang.HighOrderFTS,
             pwfts.ProbabilisticWeightedFTS]
 
 
+def get_point_multivariate_methods():
+    """Return all multivariate FTS methods por point forecasting"""
+    return [mvfts.MVFTS, wmvfts.WeightedMVFTS, cmvfts.ClusteredMVFTS]
+
+
 def get_benchmark_interval_methods():
     """Return all non FTS methods for point_to_interval forecasting"""
     return [ arima.ARIMA, quantreg.QuantileRegression]
diff --git a/pyFTS/models/hofts.py b/pyFTS/models/hofts.py
index b464eeb..917f32e 100644
--- a/pyFTS/models/hofts.py
+++ b/pyFTS/models/hofts.py
@@ -36,6 +36,48 @@ class HighOrderFLRG(flrg.FLRG):
         return len(self.RHS)
 
 
+class WeightedHighOrderFLRG(flrg.FLRG):
+    """Weighted High Order Fuzzy Logical Relationship Group"""
+
+    def __init__(self, order, **kwargs):
+        super(WeightedHighOrderFLRG, self).__init__(order, **kwargs)
+        self.LHS = []
+        self.RHS = {}
+        self.count = 0.0
+        self.strlhs = ""
+        self.w = None
+
+    def append_rhs(self, fset, **kwargs):
+        if fset not in self.RHS:
+            self.RHS[fset] = 1.0
+        else:
+            self.RHS[fset] += 1.0
+        self.count += 1.0
+
+    def append_lhs(self, c):
+        self.LHS.append(c)
+
+    def weights(self):
+        if self.w is None:
+            self.w = np.array([self.RHS[c] / self.count for c in self.RHS.keys()])
+        return self.w
+
+    def get_midpoint(self, sets):
+        mp = np.array([sets[c].centroid for c in self.RHS.keys()])
+        return mp.dot(self.weights())
+
+    def __str__(self):
+        _str = ""
+        for k in self.RHS.keys():
+            _str += ", " if len(_str) > 0 else ""
+            _str += k + " (" + str(round(self.RHS[k] / self.count, 3)) + ")"
+
+        return self.get_key() + " -> " + _str
+
+    def __len__(self):
+        return len(self.RHS)
+
+
 class HighOrderFTS(fts.FTS):
     """Conventional High Order Fuzzy Time Series"""
     def __init__(self, **kwargs):
@@ -145,7 +187,6 @@ class HighOrderFTS(fts.FTS):
         else:
             self.generate_flrg_fuzzyfied(data)
 
-
     def forecast(self, ndata, **kwargs):
 
         explain = kwargs.get('explain', False)
@@ -181,7 +222,6 @@ class HighOrderFTS(fts.FTS):
                         if explain:
                             print("\t {} -> {} (Naïve)\t Midpoint: {}\n".format(str(flrg.LHS), flrg.LHS[-1],
                                                                                             mp))
-
                 else:
                     flrg = self.flrgs[flrg.get_key()]
                     mp = flrg.get_midpoint(self.partitioner.sets)
@@ -197,3 +237,38 @@ class HighOrderFTS(fts.FTS):
                 print("Deffuzyfied value: {} \n".format(final))
 
         return ret
+
+
+class WeightedHighOrderFTS(HighOrderFTS):
+    """Weighted High Order Fuzzy Time Series"""
+    def __init__(self, **kwargs):
+        super(WeightedHighOrderFTS, self).__init__(**kwargs)
+        self.name = "Weighted High Order FTS"
+        self.shortname = "WHOFTS"
+
+    def generate_lhs_flrg_fuzzyfied(self, sample, explain=False):
+        lags = {}
+
+        flrgs = []
+
+        for ct, o in enumerate(self.lags):
+            lags[ct] = sample[o-1]
+
+            if explain:
+                print("\t (Lag {}) {} -> {} \n".format(o, sample[o-1], lhs))
+
+        root = tree.FLRGTreeNode(None)
+
+        tree.build_tree_without_order(root, lags, 0)
+
+        # Trace the possible paths
+        for p in root.paths():
+            flrg = WeightedHighOrderFLRG(self.order)
+            path = list(reversed(list(filter(None.__ne__, p))))
+
+            for lhs in path:
+                flrg.append_lhs(lhs)
+
+            flrgs.append(flrg)
+
+        return flrgs
diff --git a/pyFTS/tests/general.py b/pyFTS/tests/general.py
index c43f393..ac078b7 100644
--- a/pyFTS/tests/general.py
+++ b/pyFTS/tests/general.py
@@ -25,7 +25,7 @@ p = Grid.GridPartitioner(data=dataset, npart=20)
 
 print(p)
 
-model = hofts.HighOrderFTS(partitioner=p, order=2)
+model = hofts.WeightedHighOrderFTS(partitioner=p, order=2)
 
 model.fit(dataset) #[22, 22, 23, 23, 24])