diff --git a/ComputingTasks.cc b/ComputingTasks.cc
index 72a6f75bfbedc70abca040530c6c0211dd9472a2..a9745cf97ad0e7497293d4b5ecfea067cf0dfce2 100644
--- a/ComputingTasks.cc
+++ b/ComputingTasks.cc
@@ -219,8 +219,20 @@ VarModelStatement::writeOutput(ostream &output, const string &basename, bool min
   options_list.writeOutput(output);
   symbol_list.writeOutput("options_.var.var_list_", output);
   output << "M_.var." << name << " = options_.var;" << endl
-         << "clear options_.var;" << endl
-         << "writeVarExpectationFunction('" << name << "');" << endl;
+         << "clear options_.var;" << endl;
+}
+
+void
+VarModelStatement::createVarModelMFunction(ostream &output, const map<string, set<int> > &var_expectation_functions_to_write) const
+{
+  if (var_expectation_functions_to_write.find(name) == var_expectation_functions_to_write.end())
+    return;
+
+  stringstream ss;
+  set<int> horizons = var_expectation_functions_to_write.find(name)->second;
+  for (set<int>::const_iterator it = horizons.begin(); it != horizons.end(); it++)
+    ss << *it << " ";
+  output << "writeVarExpectationFunction('" << name << "', [" << ss.rdbuf() << "]);" << endl;
 }
 
 StochSimulStatement::StochSimulStatement(const SymbolList &symbol_list_arg,
diff --git a/ComputingTasks.hh b/ComputingTasks.hh
index 7bda35b00bb67d7893377d988ec149012ce6b8f5..651b28ea71519f0fef3a4825b45f8f3688731ec4 100644
--- a/ComputingTasks.hh
+++ b/ComputingTasks.hh
@@ -125,6 +125,7 @@ public:
   void getVarModelNameAndVarList(map<string, SymbolList > &var_model_info);
   virtual void checkPass(ModFileStructure &mod_file_struct, WarningConsolidation &warnings);
   virtual void writeOutput(ostream &output, const string &basename, bool minimal_workspace) const;
+  void createVarModelMFunction(ostream &output, const map<string, set<int> > &var_expectation_functions_to_write) const;
 };
 
 class ForecastStatement : public Statement
diff --git a/DynamicModel.cc b/DynamicModel.cc
index 7034ff63cda1b35ea303f91eb27a6688ad959687..a01a2970c4b440ccb3e6260b913c530bfa8b64da 100644
--- a/DynamicModel.cc
+++ b/DynamicModel.cc
@@ -1564,23 +1564,28 @@ DynamicModel::writeDynamicMFile(const string &dynamic_basename) const
 }
 
 void
-DynamicModel::writeVarExpectationCalls(ostream &output) const
+DynamicModel::fillVarExpectationFunctionsToWrite()
 {
-  map<string, set<int> > var_expectation_functions_to_write;
   for (var_expectation_node_map_t::const_iterator it = var_expectation_node_map.begin();
        it != var_expectation_node_map.end(); it++)
     var_expectation_functions_to_write[it->first.first].insert(it->first.second.second);
+}
+
+map<string, set<int> >
+DynamicModel::getVarExpectationFunctionsToWrite() const
+{
+  return var_expectation_functions_to_write;
+}
 
+void
+DynamicModel::writeVarExpectationCalls(ostream &output) const
+{
   for (map<string, set<int> >::const_iterator it = var_expectation_functions_to_write.begin();
        it != var_expectation_functions_to_write.end(); it++)
     {
       int i = 0;
-      stringstream ss;
-      for (set<int>::const_iterator it1 = it->second.begin(); it1 != it->second.end(); it1++)
-        ss << *it1 << " ";
-
       output << "dynamic_var_forecast_" << it->first << " = "
-             << "var_forecast_" << it->first << "(y, [" << ss.rdbuf() << "]);" << endl;
+             << "var_forecast_" << it->first << "(y);" << endl;
 
       for (set<int>::const_iterator it1 = it->second.begin(); it1 != it->second.end(); it1++)
         output << "dynamic_var_forecast_" << it->first << "_" << *it1 << " = "
diff --git a/DynamicModel.hh b/DynamicModel.hh
index 64c368cef9acee5f440b309ed35c1d406f546d09..ff9ddfd1a676a4d72a6e75b37943b8b1fc3dc909 100644
--- a/DynamicModel.hh
+++ b/DynamicModel.hh
@@ -209,6 +209,9 @@ private:
   //! List for each equation its block number
   vector<int> equation_block;
 
+  //! Used for var_expectation and var_model
+  map<string, set<int> > var_expectation_functions_to_write;
+
   //!Maximum lead and lag for each block on endogenous of the block, endogenous of the previous blocks, exogenous and deterministic exogenous
   vector<pair<int, int> > endo_max_leadlag_block, other_endo_max_leadlag_block, exo_max_leadlag_block, exo_det_max_leadlag_block, max_leadlag_block;
 
@@ -325,6 +328,12 @@ public:
   //! Transforms the model by removing trends specified by the user
   void detrendEquations();
 
+  //! Fill var_expectation_functions_to_write
+  void fillVarExpectationFunctionsToWrite();
+
+  //! Get var_expectation_functions_to_write
+  map<string, set<int> > getVarExpectationFunctionsToWrite() const;
+
   //! Transforms the model by replacing trend variables with a 1
   void removeTrendVariableFromEquations();
 
diff --git a/ModFile.cc b/ModFile.cc
index a6e180e45432f3210318fcb6a60925db4c5b7f71..310e6ce5285d5d3dd9ce40bd0c7e3bfc02c68c81 100644
--- a/ModFile.cc
+++ b/ModFile.cc
@@ -422,6 +422,7 @@ ModFile::transformPass(bool nostrict, bool compute_xrefs)
     }
   if (!var_model_info.empty())
     dynamic_model.setVarExpectationIndices(var_model_info);
+  dynamic_model.fillVarExpectationFunctionsToWrite();
 
   // Freeze the symbol table
   symbol_table.freeze();
@@ -831,10 +832,6 @@ ModFile::writeOutputFiles(const string &basename, bool clear_all, bool clear_glo
   for (vector<Statement *>::const_iterator it = statements.begin();
        it != statements.end(); it++)
     {
-      // Don't print var_model statements again as they were already printed
-      //      if (dynamic_cast<VarModelStatement *>(*it) != NULL)
-      //        continue;
-
       (*it)->writeOutput(mOutputFile, basename, minimal_workspace);
 
       /* Special treatment for initval block: insert initial values for the
@@ -855,6 +852,11 @@ ModFile::writeOutputFiles(const string &basename, bool clear_all, bool clear_glo
       LoadParamsAndSteadyStateStatement *lpass = dynamic_cast<LoadParamsAndSteadyStateStatement *>(*it);
       if (lpass && !no_static)
         static_model.writeAuxVarInitval(mOutputFile, oMatlabOutsideModel);
+
+      // Special treatement for Var Models
+      VarModelStatement *vms = dynamic_cast<VarModelStatement *>(*it);
+      if (vms != NULL)
+        vms->createVarModelMFunction(mOutputFile, dynamic_model.getVarExpectationFunctionsToWrite());
     }
 
   // Remove path for block option with M-files