diff --git a/src/DynareBison.yy b/src/DynareBison.yy
index 29913ab7c861ffb2889b95889876f690fc4927ce..d82aafa20ac5d61833f1329f4bb498493988686c 100644
--- a/src/DynareBison.yy
+++ b/src/DynareBison.yy
@@ -116,7 +116,7 @@ class ParsingDriver;
 %token <string> QUOTED_STRING
 %token QZ_CRITERIUM QZ_ZERO_THRESHOLD DSGE_VAR DSGE_VARLAG DSGE_PRIOR_WEIGHT TRUNCATE PIPE_E PIPE_X PIPE_P
 %token RELATIVE_IRF REPLIC SIMUL_REPLIC RPLOT SAVE_PARAMS_AND_STEADY_STATE PARAMETER_UNCERTAINTY TARGETS
-%token SHOCKS HETEROSKEDASTIC_SHOCKS SHOCK_DECOMPOSITION SHOCK_GROUPS USE_SHOCK_GROUPS SIGMA_E SIMUL SIMUL_ALGO SIMUL_SEED ENDOGENOUS_TERMINAL_PERIOD
+%token SHOCKS HETEROSKEDASTIC_SHOCKS SHOCK_DECOMPOSITION SHOCK_GROUPS USE_SHOCK_GROUPS SIMUL SIMUL_ALGO SIMUL_SEED ENDOGENOUS_TERMINAL_PERIOD
 %token SMOOTHER SMOOTHER2HISTVAL SQUARE_ROOT_SOLVER STACK_SOLVE_ALGO STEADY_STATE_MODEL SOLVE_ALGO SOLVER_PERIODS ROBUST_LIN_SOLVE
 %token STDERR STEADY STOCH_SIMUL SYLVESTER SYLVESTER_FIXED_POINT_TOL REGIMES REGIME REALTIME_SHOCK_DECOMPOSITION CONDITIONAL UNCONDITIONAL
 %token TEX RAMSEY_MODEL RAMSEY_POLICY RAMSEY_CONSTRAINTS PLANNER_DISCOUNT PLANNER_DISCOUNT_LATEX_NAME
@@ -245,7 +245,6 @@ statement : parameters
           | shocks
           | mshocks
           | heteroskedastic_shocks
-          | sigma_e
           | steady
           | check
           | simul
@@ -1322,8 +1321,6 @@ period_list : period_list COMMA INT_NUMBER
               }
             ;
 
-sigma_e : SIGMA_E EQUAL '[' triangular_matrix ']' ';' { driver.do_sigma_e(); };
-
 value_list : value_list COMMA '(' expression ')'
              {
                $$ = $1;
@@ -1358,26 +1355,6 @@ value_list : value_list COMMA '(' expression ')'
              }
            ;
 
-triangular_matrix : triangular_matrix ';' triangular_row
-                    { driver.end_of_row(); }
-                  | triangular_row
-                    { driver.end_of_row(); }
-                  ;
-
-triangular_row : triangular_row COMMA '(' expression ')'
-                 { driver.add_to_row($4); }
-               | triangular_row COMMA signed_number
-                 { driver.add_to_row_const($3); }
-               | triangular_row '(' expression ')'
-                 { driver.add_to_row($3); }
-               | triangular_row signed_number
-                 { driver.add_to_row_const($2); }
-               | '(' expression ')'
-                 { driver.add_to_row($2); }
-               | signed_number
-                 { driver.add_to_row_const($1); }
-               ;
-
 steady : STEADY ';'
          { driver.steady(); }
        | STEADY '(' steady_options_list ')' ';'
diff --git a/src/DynareFlex.ll b/src/DynareFlex.ll
index ed3825268d75fc0bc7cf6cb644df246aa1c0e688..40884519061b933f58c34067e0612ff54e97932c 100644
--- a/src/DynareFlex.ll
+++ b/src/DynareFlex.ll
@@ -41,10 +41,6 @@ using token = Dynare::parser::token;
 #define yyterminate() return Dynare::parser::token_type (0);
 
 int comment_caller, line_caller;
-/* Particular value : when sigma_e command is found
- this flag is set to 1, when command finished it is set to 0
- */
-int sigma_e = 0;
 string eofbuff;
 %}
 
@@ -142,7 +138,6 @@ DATE -?[0-9]+([ya]|m([1-9]|1[0-2])|q[1-4])
 <INITIAL>var_expectation_model {BEGIN DYNARE_STATEMENT; return token::VAR_EXPECTATION_MODEL;}
 <INITIAL>pac_model {BEGIN DYNARE_STATEMENT; return token::PAC_MODEL;}
 <INITIAL>dsample {BEGIN DYNARE_STATEMENT; return token::DSAMPLE;}
-<INITIAL>Sigma_e {BEGIN DYNARE_STATEMENT; sigma_e = 1; return token::SIGMA_E;}
 <INITIAL>planner_objective {BEGIN DYNARE_STATEMENT; return token::PLANNER_OBJECTIVE;}
 <INITIAL>ramsey_model {BEGIN DYNARE_STATEMENT; return token::RAMSEY_MODEL;}
 <INITIAL>ramsey_policy {BEGIN DYNARE_STATEMENT; return token::RAMSEY_POLICY;}
@@ -196,11 +191,7 @@ DATE -?[0-9]+([ya]|m([1-9]|1[0-2])|q[1-4])
 <INITIAL>var_remove {BEGIN DYNARE_STATEMENT; return token::VAR_REMOVE;}
 <INITIAL>resid {BEGIN DYNARE_STATEMENT; return token::RESID;}
 
-<DYNARE_STATEMENT>; {
-  if (!sigma_e)
-    BEGIN INITIAL;
-  return Dynare::parser::token_type (yytext[0]);
-}
+<DYNARE_STATEMENT>; {BEGIN INITIAL; return Dynare::parser::token_type (yytext[0]);}
 
 
  /* Begin of a Dynare block */
@@ -942,11 +933,7 @@ DATE -?[0-9]+([ya]|m([1-9]|1[0-2])|q[1-4])
 <DYNARE_STATEMENT,DYNARE_BLOCK>: {return Dynare::parser::token_type (yytext[0]);}
 <DYNARE_STATEMENT,DYNARE_BLOCK>[\(\)] {return Dynare::parser::token_type (yytext[0]);}
 <DYNARE_STATEMENT,DYNARE_BLOCK>\[ {return Dynare::parser::token_type (yytext[0]);}
-<DYNARE_STATEMENT,DYNARE_BLOCK>\] {
-  if (sigma_e)
-    sigma_e = 0;
-  return Dynare::parser::token_type (yytext[0]);
-}
+<DYNARE_STATEMENT,DYNARE_BLOCK>\] {return Dynare::parser::token_type (yytext[0]);}
 <DYNARE_STATEMENT,DYNARE_BLOCK>\+ {return token::PLUS;}
 <DYNARE_STATEMENT,DYNARE_BLOCK>-  {return token::MINUS;}
 <DYNARE_STATEMENT,DYNARE_BLOCK>\* {return token::TIMES;}
diff --git a/src/Makefile.am b/src/Makefile.am
index 36c5ff56db6cc8752d89753f5375e46e71bdf864..a26e7a951e0ec041aecf3abd7bf00f3ebc4de4b5 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -24,8 +24,6 @@ dynare_preprocessor_SOURCES = \
 	NumericalInitialization.hh \
 	Shocks.cc \
 	Shocks.hh \
-	SigmaeInitialization.cc \
-	SigmaeInitialization.hh \
 	SymbolTable.cc \
 	SymbolTable.hh \
 	SymbolList.cc \
diff --git a/src/ParsingDriver.cc b/src/ParsingDriver.cc
index 2d00e83cb0163cd1b70e11c2c69ed4bf325e1b94..320de782c6ae6b3c81b660461d66eeafcba61914 100644
--- a/src/ParsingDriver.cc
+++ b/src/ParsingDriver.cc
@@ -1322,43 +1322,6 @@ ParsingDriver::add_svar_global_identification_check()
   mod_file->addStatement(make_unique<SvarGlobalIdentificationCheckStatement>());
 }
 
-void
-ParsingDriver::do_sigma_e()
-{
-  warning("Sigma_e: this command is now deprecated and may be removed in a future version of Dynare. Please use the ''shocks'' command instead.");
-
-  try
-    {
-      mod_file->addStatement(make_unique<SigmaeStatement>(move(sigmae_matrix)));
-    }
-  catch (SigmaeStatement::MatrixFormException &e)
-    {
-      error("Sigma_e: matrix is neither upper triangular nor lower triangular");
-    }
-  sigmae_matrix.clear();
-}
-
-void
-ParsingDriver::end_of_row()
-{
-  sigmae_matrix.push_back(sigmae_row);
-  sigmae_row.clear();
-}
-
-void
-ParsingDriver::add_to_row_const(const string &v)
-{
-  sigmae_row.push_back(v.at(0) == '-' ?
-                       data_tree->AddUMinus(data_tree->AddNonNegativeConstant(v.substr(1))) :
-                       data_tree->AddNonNegativeConstant(v));
-}
-
-void
-ParsingDriver::add_to_row(expr_t v)
-{
-  sigmae_row.push_back(v);
-}
-
 void
 ParsingDriver::steady()
 {
diff --git a/src/ParsingDriver.hh b/src/ParsingDriver.hh
index 896c15f137361ae5e194af4cc9076de541474325..222496c9af192c3c04b2613325ff6138b6e5b3c1 100644
--- a/src/ParsingDriver.hh
+++ b/src/ParsingDriver.hh
@@ -40,7 +40,6 @@ class ParsingDriver;
 
 #include "ComputingTasks.hh"
 #include "Shocks.hh"
-#include "SigmaeInitialization.hh"
 #include "NumericalInitialization.hh"
 #include "DynamicModel.hh"
 
@@ -165,10 +164,6 @@ private:
   ShocksStatement::covar_and_corr_shocks_t corr_shocks;
   //! Temporary storage for values and scales of heteroskedastic_shocks
   HeteroskedasticShocksStatement::heteroskedastic_shocks_t heteroskedastic_shocks_values, heteroskedastic_shocks_scales;
-  //! Temporary storage for Sigma_e rows
-  SigmaeStatement::row_t sigmae_row;
-  //! Temporary storage for Sigma_e matrix
-  SigmaeStatement::matrix_t sigmae_matrix;
   //! Temporary storage for initval blocks
   InitOrEndValStatement::init_values_t init_values;
   /* Temporary storage for endval blocks. Uses a type that encompasses both
@@ -486,14 +481,6 @@ public:
   //! Adds a deterministic shock value
   /*! \param v a string containing a (possibly negative) numeric constant */
   void add_value(const string &v);
-  //! Writes a Sigma_e block
-  void do_sigma_e();
-  //! Ends row of Sigma_e block
-  void end_of_row();
-  //! Adds a constant element to current row of Sigma_e
-  void add_to_row_const(const string &v);
-  //! Adds an expression element to current row of Sigma_e
-  void add_to_row(expr_t v);
   //! Write a steady command
   void steady();
   //! Sets an option to a numerical value
diff --git a/src/SigmaeInitialization.cc b/src/SigmaeInitialization.cc
deleted file mode 100644
index 07cded9ba82eb16d02b5358d9384428b06487be7..0000000000000000000000000000000000000000
--- a/src/SigmaeInitialization.cc
+++ /dev/null
@@ -1,123 +0,0 @@
-/*
- * Copyright © 2003-2022 Dynare Team
- *
- * This file is part of Dynare.
- *
- * Dynare is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * Dynare is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
- */
-
-#include <utility>
-
-#include "SigmaeInitialization.hh"
-
-SigmaeStatement::SigmaeStatement(matrix_t matrix_arg) noexcept(false) :
-  matrix{move(matrix_arg)},
-  matrix_form{determineMatrixForm(matrix)}
-{
-}
-
-SigmaeStatement::MatrixForm
-SigmaeStatement::determineMatrixForm(const matrix_t &matrix) noexcept(false)
-{
-  size_t nbe;
-  int inc;
-  MatrixForm type;
-  // Checking if first or last row has one element.
-  if (matrix.front().size() == 1)
-    {
-      inc = 1;
-      nbe = 2;
-      type = MatrixForm::lower;
-    }
-  else if (matrix.back().size() == 1)
-    {
-      inc = -1;
-      nbe = matrix.front().size()-1;
-      type = MatrixForm::upper;
-    }
-  else
-    throw MatrixFormException();
-
-  // Checking if matrix is triangular (upper or lower):
-  // each row has one element more or less than the previous one
-  // and first or last one has one element.
-  for (auto ir = ++matrix.begin(); ir != matrix.end(); ++ir, nbe += inc)
-    if (ir->size() != nbe)
-      throw MatrixFormException();
-
-  return type;
-}
-
-void
-SigmaeStatement::writeOutput(ostream &output, [[maybe_unused]] const string &basename,
-                             [[maybe_unused]] bool minimal_workspace) const
-{
-  output << "M_.Sigma_e = [..." << endl;
-  for (size_t ir = 0; ir < matrix.size(); ir++)
-    {
-      for (size_t ic = 0; ic < matrix.size(); ic++)
-        {
-          size_t ic1, ir1;
-          if (ic >= ir && matrix_form == MatrixForm::upper)
-            {
-              ic1 = ic-ir;
-              ir1 = ir;
-            }
-          else if (ic < ir && matrix_form == MatrixForm::upper)
-            {
-              ic1 = ir-ic;
-              ir1 = ic;
-            }
-          else if (ic > ir && matrix_form == MatrixForm::lower)
-            {
-              ic1 = ir;
-              ir1 = ic;
-            }
-          else // ic <= ir && matrix_form == MatrixForm::lower
-            {
-              ic1 = ic;
-              ir1 = ir;
-            }
-
-          matrix[ir1][ic1]->writeOutput(output);
-          output << " ";
-        }
-      output << ";..." << endl;
-    }
-  output << "];" << endl;
-}
-
-void
-SigmaeStatement::writeJsonOutput(ostream &output) const
-{
-  output << R"({"statementName": "Sigma_e", "value": [)";
-  for (bool printed_something{false};
-       const auto &row : matrix)
-    {
-      if (exchange(printed_something, true))
-        output << ", ";
-      output << "[";
-      for (bool printed_something2{false};
-           auto elem : row)
-        {
-          if (exchange(printed_something2, true))
-            output << ", ";
-          output << '"';
-          elem->writeJsonOutput(output, {}, {});
-          output << '"';
-        }
-      output << "]";
-    }
-  output << "] }";
-}
diff --git a/src/SigmaeInitialization.hh b/src/SigmaeInitialization.hh
deleted file mode 100644
index 1d5bf1c2140d1ada3317bc92601769efe9d265d5..0000000000000000000000000000000000000000
--- a/src/SigmaeInitialization.hh
+++ /dev/null
@@ -1,64 +0,0 @@
-/*
- * Copyright © 2003-2020 Dynare Team
- *
- * This file is part of Dynare.
- *
- * Dynare is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * Dynare is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
- */
-
-#ifndef _SIGMAEINITIALIZATION_HH
-#define _SIGMAEINITIALIZATION_HH
-
-#include <string>
-#include <vector>
-
-#include "ExprNode.hh"
-#include "Statement.hh"
-
-//! Stores a Sigma_e statement
-class SigmaeStatement : public Statement
-{
-public:
-  //! Matrix form (lower or upper triangular) enum
-  enum class MatrixForm
-    {
-     lower, //!< Lower triangular matrix
-     upper //!< Upper triangular matrix
-    };
-  //! Type of a matrix row
-  using row_t = vector<expr_t>;
-  //! Type of a complete matrix
-  using matrix_t = vector<row_t>;
-
-  //! An exception indicating that a matrix is neither upper triangular nor lower triangular
-  class MatrixFormException
-  {
-  };
-private:
-  //! The matrix
-  const matrix_t matrix;
-  //! Matrix form (lower or upper)
-  const MatrixForm matrix_form;
-
-  //! Returns the type (upper or lower triangular) of a given matrix
-  /*! Throws an exception if it is neither upper triangular nor lower triangular */
-  static MatrixForm determineMatrixForm(const matrix_t &matrix) noexcept(false);
-
-public:
-  explicit SigmaeStatement(matrix_t matrix_arg) noexcept(false);
-  void writeOutput(ostream &output, const string &basename, bool minimal_workspace) const override;
-  void writeJsonOutput(ostream &output) const override;
-};
-
-#endif