diff --git a/preprocessor/DynamicModel.cc b/preprocessor/DynamicModel.cc
index 51ca2bebb..a4a6ad799 100644
--- a/preprocessor/DynamicModel.cc
+++ b/preprocessor/DynamicModel.cc
@@ -2678,6 +2678,14 @@ DynamicModel::computeParamsDerivatives()
 
       int param = it->second;
 
+      for (int eq = 0; eq < (int) equations.size(); eq++)
+        {
+          NodeID d1 = equations[eq]->getDerivative(param);
+          if (d1 == Zero)
+            continue;
+          residuals_params_derivatives[make_pair(eq, param)] = d1;
+        }
+
       for (first_derivatives_type::const_iterator it2 = first_derivatives.begin();
            it2 != first_derivatives.end(); it2++)
         {
@@ -2688,7 +2696,7 @@ DynamicModel::computeParamsDerivatives()
           NodeID d2 = d1->getDerivative(param);
           if (d2 == Zero)
             continue;
-          params_derivatives[make_pair(eq, make_pair(var, param))] = d2;
+          jacobian_params_derivatives[make_pair(eq, make_pair(var, param))] = d2;
         }
     }
 }
@@ -2699,15 +2707,20 @@ DynamicModel::computeParamsDerivativesTemporaryTerms()
   map<NodeID, int> reference_count;
   params_derivs_temporary_terms.clear();
 
-  for (second_derivatives_type::iterator it = params_derivatives.begin();
-       it != params_derivatives.end(); it++)
+  for (first_derivatives_type::iterator it = residuals_params_derivatives.begin();
+       it != residuals_params_derivatives.end(); it++)
+    it->second->computeTemporaryTerms(reference_count, params_derivs_temporary_terms, true);
+
+  for (second_derivatives_type::iterator it = jacobian_params_derivatives.begin();
+       it != jacobian_params_derivatives.end(); it++)
     it->second->computeTemporaryTerms(reference_count, params_derivs_temporary_terms, true);
 }
 
 void
 DynamicModel::writeParamsDerivativesFile(const string &basename) const
   {
-    if (!params_derivatives.size())
+    if (!residuals_params_derivatives.size()
+        && !jacobian_params_derivatives.size())
       return;
 
     string filename = basename + "_params_derivs.m";
@@ -2719,7 +2732,7 @@ DynamicModel::writeParamsDerivativesFile(const string &basename) const
         cerr << "ERROR: Can't open file " << filename << " for writing" << endl;
         exit(EXIT_FAILURE);
       }
-    paramsDerivsFile << "function gp = " << basename << "_params_derivs(y, x, params, it_)" << endl
+    paramsDerivsFile << "function [rp, gp] = " << basename << "_params_derivs(y, x, params, it_)" << endl
     << "%" << endl
     << "% Warning : this file is generated automatically by Dynare" << endl
     << "%           from model file (.mod)" << endl << endl;
@@ -2727,11 +2740,30 @@ DynamicModel::writeParamsDerivativesFile(const string &basename) const
 
     writeTemporaryTerms(params_derivs_temporary_terms, paramsDerivsFile, oMatlabDynamicModel);
 
+    // Write parameter derivative
+    paramsDerivsFile << "rp = zeros(" << equation_number() << ", "
+                     << symbol_table.param_nbr() << ");" << endl;
+
+    for (first_derivatives_type::const_iterator it = residuals_params_derivatives.begin();
+         it != residuals_params_derivatives.end(); it++)
+      {
+        int eq = it->first.first;
+        int param = it->first.second;
+        NodeID d1 = it->second;
+
+        int param_col = symbol_table.getTypeSpecificID(getSymbIDByDerivID(param)) + 1;
+
+        paramsDerivsFile << "rp(" << eq+1 << ", " << param_col << ") = ";
+        d1->writeOutput(paramsDerivsFile, oMatlabDynamicModel, params_derivs_temporary_terms);
+        paramsDerivsFile << ";" << endl;
+      }
+
+    // Write jacobian derivatives
     paramsDerivsFile << "gp = zeros(" << equation_number() << ", " << dynJacobianColsNbr << ", "
     << symbol_table.param_nbr() << ");" << endl;
 
-    for (second_derivatives_type::const_iterator it = params_derivatives.begin();
-         it != params_derivatives.end(); it++)
+    for (second_derivatives_type::const_iterator it = jacobian_params_derivatives.begin();
+         it != jacobian_params_derivatives.end(); it++)
       {
         int eq = it->first.first;
         int var = it->first.second.first;
diff --git a/preprocessor/DynamicModel.hh b/preprocessor/DynamicModel.hh
index f1b5dbfc5..7598faa53 100644
--- a/preprocessor/DynamicModel.hh
+++ b/preprocessor/DynamicModel.hh
@@ -60,12 +60,19 @@ private:
   /*! Set by computeDerivID() and computeDynJacobianCols() */
   int dynJacobianColsNbr;
 
+  //! Derivatives of the residuals w.r. to parameters
+  /*! First index is equation number, second is parameter.
+    Only non-null derivatives are stored in the map.
+    Parameter indices are those of the getDerivID() method.
+  */
+  first_derivatives_type residuals_params_derivatives;
+
   //! Derivatives of the jacobian w.r. to parameters
   /*! First index is equation number, second is endo/exo/exo_det variable, and third is parameter.
     Only non-null derivatives are stored in the map.
     Variable and parameter indices are those of the getDerivID() method.
   */
-  second_derivatives_type params_derivatives;
+  second_derivatives_type jacobian_params_derivatives;
 
   //! Temporary terms for the file containing parameters dervicatives
   temporary_terms_type params_derivs_temporary_terms;