diff --git a/doc/manual/source/time-series.rst b/doc/manual/source/time-series.rst
index 3c5b682d5b53cc5e8970c227a4ac96bd4818b392..f7dd8d484cf447168649418e3367cd8001c9fb94 100644
--- a/doc/manual/source/time-series.rst
+++ b/doc/manual/source/time-series.rst
@@ -1013,26 +1013,34 @@ The dseries class
         input. Valid file types are ``.m``, ``.mat``, ``.csv`` and
         ``.xls/.xlsx`` (Octave only supports ``.xlsx`` files and the
         `io <https://octave.sourceforge.io/io/>`__ package from
-        Octave-Forge must be installed). A typical ``.m`` file will
-        have the following form::
+        Octave-Forge must be installed). The extension of the file
+        should be explicitly provided. A typical ``.m`` file will have
+        the following form::
 
+            FREQ__ = 4;
             INIT__ = '1994Q3';
             NAMES__ = {'azert';'yuiop'};
             TEX__ = {'azert';'yuiop'};
+            TAGS__ = struct()
+            DATA__ = {}
 
             azert = randn(100,1);
             yuiop = randn(100,1);
 
         If a ``.mat`` file is used instead, it should provide the same
-        informations. Note that the ``INIT__`` variable can be either
-        a ``dates`` object or a string which could be used to
-        instantiate the same ``dates`` object. If ``INIT__`` is not
-        provided in the ``.mat`` or ``.m`` file, the initial is by
-        default set equal to ``dates('1Y')``. If a second input
-        argument is passed to the constructor, ``dates`` object
-        *INITIAL_DATE*, the initial date defined in *FILENAME* is
-        reset to *INITIAL_DATE*. This is typically usefull if
-        ``INIT__`` is not provided in the data file.
+        informations, except that the data should not be given as a
+        set of vectors, but as a single matrix of doubles named
+        ``DATA__``. This array should have as many columns as elements
+        in ``NAMES__`` (the number of variables). Note that the
+        ``INIT__`` variable can be either a ``dates`` object or a
+        string which could be used to instantiate the same ``dates``
+        object. If ``INIT__`` is not provided in the ``.mat`` or
+        ``.m`` file, the initial is by default set equal to
+        ``dates('1Y')``. If a second input argument is passed to the
+        constructor, ``dates`` object *INITIAL_DATE*, the initial date
+        defined in *FILENAME* is reset to *INITIAL_DATE*. This is
+        typically usefull if ``INIT__`` is not provided in the data
+        file.
 
     .. construct:: dseries (DATA_MATRIX[,INITIAL_DATE[,LIST_OF_NAMES[,TEX_NAMES]]])
                    dseries (DATA_MATRIX[,RANGE_OF_DATES[,LIST_OF_NAMES[,TEX_NAMES]]])
@@ -1062,20 +1070,20 @@ The dseries class
 
        Creates a ``dseries`` object given the MATLAB Table provided as the sole
        argument. It is assumed that the first column of the table contains the
-       dates of the ``dseries`` and the first row contains the names. NB: This
+       dates of the ``dseries`` and the first row contains the names. This
        feature is not available under Octave or MATLAB R2013a or earlier.
 
-    *Example*
+       *Example*
 
-        Various ways to create a ``dseries`` object::
+       Various ways to create a ``dseries`` object::
 
-            do1 = dseries(1999Q3);
-            do2 = dseries('filename.csv');
-            do3 = dseries([1; 2; 3], 1999Q3, {'var123'}, {'var_{123}'});
+         do1 = dseries(1999Q3);
+         do2 = dseries('filename.csv');
+         do3 = dseries([1; 2; 3], 1999Q3, {'var123'}, {'var_{123}'});
 
-            >> do1 = dseries(dates('1999Q3'));
-            >> do2 = dseries('filename.csv');
-            >> do3 = dseries([1; 2; 3], dates('1999Q3'), {'var123'}, {'var_{123}'});
+         >> do1 = dseries(dates('1999Q3'));
+         >> do2 = dseries('filename.csv');
+         >> do3 = dseries([1; 2; 3], dates('1999Q3'), {'var123'}, {'var_{123}'});
 
 
     One can easily create subsamples from a ``dseries`` object using
@@ -1087,10 +1095,13 @@ The dseries class
     to the last observation, then ``ds(d)`` instantiates a new
     ``dseries`` object containing the subsample defined by ``d``.
 
-    A list of the available methods, by alphabetical order, is given below.
+    A list of the available methods, by alphabetical order, is given
+    below. As in the previous section the in place modifications
+    versions of the methods are postfixed with an underscore.
 
 
-    .. dseriesmethod:: A = abs(B)
+    .. dseriesmethod:: A = abs (B)
+                       abs_ (B)
 
         |br| Overloads the ``abs()`` function for ``dseries``
         objects. Returns the absolute value of the variables in
@@ -1121,7 +1132,8 @@ The dseries class
                 1973Q3 | 0.99791 | 0.22677
 
 
-    .. dseriesmethod:: [A, B] = align(A, B)
+    .. dseriesmethod:: [A, B] = align (A, B)
+                       align_ (A, B)
 
         If ``dseries`` objects ``A`` and ``B`` are defined on
         different time ranges, this function extends ``A`` and/or
@@ -1160,8 +1172,33 @@ The dseries class
                 2001Q1 | 0.17813
                 2001Q2 | 0.12801
 
+                >> ts0 = dseries(rand(5,1),dates('2000Q1')); % 2000Q1 -> 2001Q1
+                >> ts1 = dseries(rand(3,1),dates('2000Q4')); % 2000Q4 -> 2001Q2
+                >> align_(ts0, ts1);                         % 2000Q1 -> 2001Q2
+                >> ts1
+
+                ts1 is a dseries object:
 
-    .. dseriesmethod:: B = baxter_king_filter(A, hf, lf, K)
+                       | Variable_1
+                2000Q1 | NaN
+                2000Q2 | NaN
+                2000Q3 | NaN
+                2000Q4 | 0.66653
+                2001Q1 | 0.17813
+                2001Q2 | 0.12801
+
+
+    .. dseriesmethod:: C = backcast (A, B[, diff])
+                       backcast_ (A, B[, diff])
+
+        Backcasts ``dseries`` object ``A`` with ``dseries`` object B's
+        growth rates (except if the last optional argument, ``diff``,
+        is true in which case first differences are used). Both
+        ``dseries`` objects must have the same frequency.
+
+
+    .. dseriesmethod:: B = baxter_king_filter (A, hf, lf, K)
+                       baxter_king_filter_ (A, hf, lf, K)
 
         |br| Implementation of the *Baxter and King* (1999) band pass
         filter for ``dseries`` objects. This filter isolates business
@@ -1206,7 +1243,17 @@ The dseries class
                 set(gca,'XTickLabel',strings(ts1.dates(id)));
 
 
-    .. dseriesmethod:: C = chain(A, B)
+    .. dseriesmethod:: B = center (A[, geometric])
+                       center_ (A[, geometric])
+
+       |br| Centers variables in ``dseries`` object ``A`` around their
+       arithmetic means, except if the optional argument ``geometric``
+       is set equal to ``true`` in which case all the variables are
+       divided by their geometric means.
+
+
+    .. dseriesmethod:: C = chain (A, B)
+                       chain_ (A, B)
 
         |br| Merge two ``dseries`` objects along the time
         dimension. The two objects must have the same number of
@@ -1252,7 +1299,7 @@ The dseries class
                 1951Q2 | 6
 
 
-    .. dseriesmethod:: [error_flag, message ] = check(A)
+    .. dseriesmethod:: [error_flag, message ] = check (A)
 
         |br| Sanity check of ``dseries`` object ``A``. Returns ``1``
         if there is an error, ``0`` otherwise. The second output
@@ -1260,7 +1307,68 @@ The dseries class
         error.
 
 
-    .. dseriesmethod:: B = cumprod(A[, d[, v]])
+    .. dseriesmethod:: B = copy (A)
+
+       |br| Returns a copy of ``A``. If an inplace modification method
+       is applied to ``A``, object ``B`` will not be affected. Note
+       that if ``A`` is assigned to ``C``, ``C = A``, then any in
+       place modification method applied to ``A`` will change ``C``.
+
+       *Example*
+
+            ::
+
+               >> a = dseries(randn(5,1))
+
+               a is a dseries object:
+
+                  | Variable_1
+               1Y | -0.16936
+               2Y | -1.1451
+               3Y | -0.034331
+               4Y | -0.089042
+               5Y | -0.66997
+
+               >> b = copy(a);
+               >> c = a;
+               >> a.abs();
+               >> a.abs_();
+               >> a
+
+               a is a dseries object:
+
+                  | Variable_1
+               1Y | 0.16936
+               2Y | 1.1451
+               3Y | 0.034331
+               4Y | 0.089042
+               5Y | 0.66997
+
+               >> b
+
+               b is a dseries object:
+
+                  | Variable_1
+               1Y | -0.16936
+               2Y | -1.1451
+               3Y | -0.034331
+               4Y | -0.089042
+               5Y | -0.66997
+
+               >> c
+
+               c is a dseries object:
+
+                  | Variable_1
+               1Y | 0.16936
+               2Y | 1.1451
+               3Y | 0.034331
+               4Y | 0.089042
+               5Y | 0.66997
+
+
+   .. dseriesmethod:: B = cumprod (A[, d[, v]])
+                      cumprod_ (A[, d[, v]])
 
         |br| Overloads the MATLAB/Octave ``cumprod`` function for
         ``dseries`` objects. The cumulated product cannot be computed
@@ -1322,7 +1430,8 @@ The dseries class
                 7Y | 50.2655
 
 
-    .. dseriesmethod:: B = cumsum(A[, d[, v]])
+    .. dseriesmethod:: B = cumsum (A[, d[, v]])
+                       cumsum (A[, d[, v]])
 
         |br| Overloads the MATLAB/Octave ``cumsum`` function for
         ``dseries`` objects. The cumulated sum cannot be computed if
@@ -1393,14 +1502,46 @@ The dseries class
                 10Y | 10.1416
 
 
-    .. dseriesmethod:: C = eq(A, B)
+    .. dseriesmethod:: B = detrend (A, m)
+                       dentrend_ (A, m)
+
+        |br| Detrends ``dseries`` object ``A`` with a fitted
+        polynomial of order ``m``. Note that each variable is
+        detrended with a different polynomial.
+
+
+    .. dseriesmethod:: B = diff (A)
+                       diff_ (A)
+
+        |br| Returns the first difference of ``dseries`` object ``A``.
+
+
+    .. datesmethod:: disp (A)
+
+        |br| Overloads the MATLAB/Octave disp function for ``dseries`` object.
+
+
+    .. datesmethod:: display (A)
+
+        |br| Overloads the MATLAB/Octave display function for
+        ``dseries`` object. ``display`` is the function called by
+        MATLAB to print the content of an object if a semicolon is
+        missing at the end of a MATLAB statement. If the ``dseries``
+        object is defined over a too large time span, only the first
+        and last periods will be printed. If the ``dseries`` object
+        contains too many variables, only the first and last variables
+        will be printed. If all the periods and variables are
+        required, the ``disp`` method should be used instead.
+
+
+    .. dseriesmethod:: C = eq (A, B)
 
         |br| Overloads the MATLAB/Octave ``eq`` (equal, ``==``)
         operator. ``dseries`` objects ``A`` and ``B`` must have the
         same number of observations (say, :math:`T`) and variables
         (:math:`N`). The returned argument is a :math:`T \times N`
-        matrix of zeros and ones. Element :math:`(i,j)` of ``C`` is
-        equal to ``1`` if and only if observation :math:`i` for
+        matrix of logicals. Element :math:`(i,j)` of ``C`` is
+        equal to ``true`` if and only if observation :math:`i` for
         variable :math:`j` in ``A`` and ``B`` are the same.
 
         *Example*
@@ -1413,28 +1554,17 @@ The dseries class
 
                 ans =
 
-                     1
-                     0
-                     1
-
-
-    .. dseriesmethod:: B = exp(A)
-
-        |br| Overloads the MATLAB/Octave ``exp`` function for
-        ``dseries`` objects.
-
-        *Example*
+                   3x1 logical array
 
-            ::
-
-                >> ts0 = dseries(rand(10,1));
-                >> ts1 = ts0.exp();
+                    1
+                    0
+                    1
 
 
-    .. dseriesmethod:: l = exist(A, varname)
+    .. dseriesmethod:: l = exist (A, varname)
 
-        |br| Tests if variable exists in ``dseries`` object ``A``. Returns
-        ``1`` (true) iff variable exists in ``A``.
+        |br| Tests if variable ``varname``  exists in ``dseries`` object ``A``. Returns
+        ``true`` iff variable exists in ``A``.
 
         *Example*
 
@@ -1445,16 +1575,34 @@ The dseries class
 
                 ans =
 
-                     1
+                   logical
+
+                    1
 
                 >> ts.exist('Variable_2')
 
                 ans =
 
-                     0
+                   logical
+
+                    0
+
+
+    .. dseriesmethod:: B = exp (A)
+                       exp_ (A)
+
+        |br| Overloads the MATLAB/Octave ``exp`` function for
+        ``dseries`` objects.
+
+        *Example*
+
+            ::
 
+                >> ts0 = dseries(rand(10,1));
+                >> ts1 = ts0.exp();
 
-    .. dseriesmethod:: C = extract(A, B[, ...])
+
+    .. dseriesmethod:: C = extract (A, B[, ...])
 
         |br| Extracts some variables from a ``dseries`` object ``A``
         and returns a ``dseries`` object ``C``. The input arguments
@@ -1472,13 +1620,15 @@ The dseries class
 
                 >> ts0 = dseries(ones(100,10));
                 >> ts1 = ts0{'Variable_1','Variable_2','Variable_3'};
-                >> ts2 = ts0{'Variable_@1,2,3@'}
-                >> ts3 = ts0{'Variable_[1-3]$'}
+                >> ts2 = ts0{'Variable_@1,2,3@'};
+                >> ts3 = ts0{'Variable_[1-3]$'};
                 >> isequal(ts1,ts2) && isequal(ts1,ts3)
 
                 ans =
 
-                     1
+                   logical
+
+                    1
 
             It is possible to use up to two implicit loops to select variables::
 
@@ -1499,7 +1649,17 @@ The dseries class
                 1973Q4 | -0.03705 | -0.35899  | 0.85838   | -1.4675  | -2.1666  | -0.62032
 
 
-    .. dseriesmethod:: f = freq(B)
+    .. dseriesmethod:: f = firstdate (A)
+
+       |br| Returns the first period in ``dseries`` object ``A``.
+
+
+    .. dseriesmethod:: f = firstobservedperiod (A)
+
+       |br| Returns the first period where all the variables in ``dseries`` object ``A`` are observed (non NaN).
+
+
+    .. dseriesmethod:: f = frequency (B)
 
         |br| Returns the frequency of the variables in ``dseries`` object ``B``.
 
@@ -1508,14 +1668,14 @@ The dseries class
             ::
 
                 >> ts = dseries(randn(3,2),'1973Q1');
-                >> ts.freq
+                >> ts.frequency
 
                 ans =
 
                      4
 
 
-    .. dseriesmethod:: D = horzcat(A, B[, ...])
+    .. dseriesmethod:: D = horzcat (A, B[, ...])
 
         |br| Overloads the ``horzcat`` MATLAB/Octave’s method for
         ``dseries`` objects. Returns a ``dseries`` object ``D``
@@ -1550,7 +1710,8 @@ The dseries class
                 1952Q1 | NaN     | NaN      | 0.50946
 
 
-    .. dseriesmethod:: B = hpcycle(A[, lambda])
+    .. dseriesmethod:: B = hpcycle (A[, lambda])
+                       hpcycle_ (A[, lambda])
 
         |br| Extracts the cycle component from a ``dseries`` ``A``
         object using the *Hodrick and Prescott (1997)* filter and
@@ -1590,7 +1751,8 @@ The dseries class
                 set(gca,'XTickLabel',strings(ts.dates(id)));
 
 
-    .. dseriesmethod:: B = hptrend(A[, lambda])
+    .. dseriesmethod:: B = hptrend (A[, lambda])
+                       hptrend_ (A[, lambda])
 
         |br| Extracts the trend component from a ``dseries`` A object
         using the *Hodrick and Prescott (1997)* filter and returns a
@@ -1618,20 +1780,7 @@ The dseries class
                 set(gca,'XTickLabel',strings(ts0.dates(id)));
 
 
-    .. dseriesmethod:: f = init(B)
-
-        |br| Returns the initial date in ``dseries`` object ``B``.
-
-        *Example*
-
-            ::
-
-                >> ts = dseries(randn(3,2),'1973Q1');
-                >> ts.init
-                ans = <dates: 1973Q1>
-
-
-    .. dseriesmethod:: C = insert(A, B, I)
+    .. dseriesmethod:: C = insert (A, B, I)
 
         |br| Inserts variables contained in ``dseries`` object ``B``
         in ``dseries`` object ``A`` at positions specified by integer
@@ -1642,7 +1791,7 @@ The dseries class
         defined over the same time ranges, but it is assumed that they
         have common frequency.
 
-        :ex:
+        *Example*
 
             ::
 
@@ -1666,29 +1815,50 @@ The dseries class
                 1950Q2 | 1   | 1     | 3.1416 | 1      | 1.7725      | 1
 
 
-    .. dseriesmethod:: B = isempty(A)
+    .. dseriesmethod:: B = isempty (A)
 
-    |br| Overloads the MATLAB/octave’s ``isempty`` function. Returns
-    ``1`` if ``dseries`` object ``A`` is empty, ``0`` otherwise.
+       |br| Overloads the MATLAB/octave’s ``isempty`` function. Returns
+       ``true`` if ``dseries`` object ``A`` is empty.
 
 
-    .. dseriesmethod:: C = isequal(A,B)
+    .. dseriesmethod:: C = isequal (A, B)
 
         |br| Overloads the MATLAB/octave’s ``isequal`` function. Returns
-        ``1`` if ``dseries`` objects ``A`` and ``B`` are identical, ``0``
-        otherwise.
+        ``true`` if ``dseries`` objects ``A`` and ``B`` are identical.
+
+
+    .. dseriesmethod:: C = isinf (A)
+
+        |br| Overloads the MATLAB/octave’s ``isinf`` function. Returns
+        a logical array, with element ``(i,j)`` equal to ``true`` if and
+        only if variable ``j`` is finite in period ``A.dates(i)``.
+
+
+    .. dseriesmethod:: C = isnan (A)
+
+        |br| Overloads the MATLAB/octave’s ``isnan`` function. Returns
+        a logical array, with element ``(i,j)`` equal to ``true`` if and
+        only if variable ``j`` isn't NaN in period ``A.dates(i)``.
 
 
-    .. dseriesmethod:: B = lag(A[, p])
+    .. dseriesmethod:: C = isreal (A)
 
-        Returns lagged time series. Default value of ``p``, the number
+        |br| Overloads the MATLAB/octave’s ``isreal`` function. Returns
+        a logical array, with element ``(i,j)`` equal to ``true`` if and
+        only if variable ``j`` is real in period ``A.dates(i)``.
+
+
+    .. dseriesmethod:: B = lag (A[, p])
+                       lag_ (A[, p])
+
+        |br| Returns lagged time series. Default value of integer scalar ``p``, the number
         of lags, is ``1``.
 
         *Example*
 
             ::
 
-                >> ts0 = dseries(transpose(1:4),'1950Q1')
+                >> ts0 = dseries(transpose(1:4), '1950Q1')
 
                 ts0 is a dseries object:
 
@@ -1702,7 +1872,7 @@ The dseries class
 
                 ts1 is a dseries object:
 
-                           | lag(Variable_1,1)
+                           | Variable_1
                     1950Q1 | NaN
                     1950Q2 | 1
                     1950Q3 | 2
@@ -1712,7 +1882,7 @@ The dseries class
 
                 ts2 is a dseries object:
 
-                       | lag(Variable_1,2)
+                       | Variable_1
                 1950Q1 | NaN
                 1950Q2 | NaN
                 1950Q3 | 1
@@ -1726,7 +1896,7 @@ The dseries class
 
                 ans is a dseries object:
 
-                       | lag(Variable_1,1)
+                       | Variable_1
                 1950Q1 | NaN
                 1950Q2 | 1
                 1950Q3 | 2
@@ -1738,32 +1908,39 @@ The dseries class
 
                 ans is a dseries object:
 
-                       | lag(Variable_1,1)
+                       | Variable_1
                 1950Q1 | NaN
                 1950Q2 | 1
                 1950Q3 | 2
                 1950Q4 | 3
 
 
-    .. dseriesmethod:: l = last(B)
+    .. dseriesmethod:: l = lastdate (B)
 
-        |br| Returns the last date in ``dseries`` object ``B``.
+        |br| Returns the last period in ``dseries`` object ``B``.
 
         *Example*
 
             ::
 
                 >> ts = dseries(randn(3,2),'1973Q1');
-                >> ts.last
+                >> ts.lastdate()
+
                 ans = <dates: 1973Q3>
 
 
-    .. dseriesmethod:: B = lead(A[, p])
+    .. dseriesmethod:: f = lastobservedperiod (A)
 
-        |br| Returns lead time series. Default value of ``p``, the
-        number of leads, is ``1``. As in the ``lag`` method, the
-        ``dseries`` class overloads the parenthesis so that
-        ``ts.lead(p)`` is equivalent to ``ts(p)``.
+       |br| Returns the last period where all the variables in ``dseries`` object ``A`` are observed (non NaN).
+
+
+    .. dseriesmethod:: B = lead (A[, p])
+                       lead_ (A[, p])
+
+        |br| Returns lead time series. Default value of integer scalar
+        ``p``, the number of leads, is ``1``. As in the ``lag``
+        method, the ``dseries`` class overloads the parenthesis so
+        that ``ts.lead(p)`` is equivalent to ``ts(p)``.
 
         *Example*
 
@@ -1774,7 +1951,7 @@ The dseries class
 
                 ts1 is a dseries object:
 
-                       | lead(Variable_1,1)
+                       | Variable_1
                 1950Q1 | 2
                 1950Q2 | 3
                 1950Q3 | 4
@@ -1784,7 +1961,7 @@ The dseries class
 
                 ts2 is a dseries object:
 
-                       | lead(Variable_1,2)
+                       | Variable_1
                 1950Q1 | 3
                 1950Q2 | 4
                 1950Q3 | NaN
@@ -1815,7 +1992,29 @@ The dseries class
         in the ``model`` block is zero, but the residuals are non
         zero).
 
-    .. dseriesmethod:: B = log(A)
+
+    .. dseriesmethod:: B = lineartrend (A)
+
+        |br| Returns a linear trend centered on 0, the length of the
+        trend is given by the size of ``dseries`` object ``A`` (the
+        number of periods).
+
+        *Example*
+
+            ::
+
+               >> ts = dseries(ones(3,1));
+               >> ts.lineartrend()
+
+               ans =
+
+                    -1
+                     0
+                     1
+
+
+    .. dseriesmethod:: B = log (A)
+                       log_ (A)
 
         |br| Overloads the MATLAB/Octave ``log`` function for
         ``dseries`` objects.
@@ -1827,8 +2026,23 @@ The dseries class
                 >> ts0 = dseries(rand(10,1));
                 >> ts1 = ts0.log();
 
+    .. dseriesmethod:: B = mdiff (A)
+                       mdiff_ (A)
+
+       |br| Computes monthly growth rates of variables in
+       ``dseries`` object ``A``.
+
 
-    .. dseriesmethod:: C = merge(A, B)
+    .. dseriesmethod:: B = mean (A[, geometric])
+
+        |br| Overloads the MATLAB/Octave ``mean`` function for
+        ``dseries`` objects. Returns the mean of each variable in
+        ``dseries`` object ``A``. If the second argument is ``true``
+        the geometric mean is computed, otherwise (default) the
+        arithmetic mean is reported.
+
+
+    .. dseriesmethod:: C = merge (A, B[, legacy])
 
         |br| Merges two ``dseries`` objects ``A`` and ``B`` in
         ``dseries`` object ``C``. Objects ``A`` and ``B`` need to have
@@ -1838,71 +2052,75 @@ The dseries class
         select the variable ``x`` as defined in the second input
         argument, ``B``, except for the NaN elements in ``B`` if
         corresponding elements in ``A`` (ie same periods) are well
-        defined numbers.
+        defined numbers. This behaviour can be changed by setting the
+        optional argument ``legacy`` equal to true, in which case the
+        second variable overwrites the first one even if the second
+        variable has NaNs.
 
         *Example*
 
             ::
 
-                >> ts0 = dseries(rand(3,2),'1950Q1',{'A1';'A2'})
+               >> ts0 = dseries(rand(3,2),'1950Q1',{'A1';'A2'})
 
-                ts0 is a dseries object:
+               ts0 is a dseries object:
 
-                       | A1       | A2
-                1950Q1 | 0.42448  | 0.92477
-                1950Q2 | 0.60726  | 0.64208
-                1950Q3 | 0.070764 | 0.1045
+                      | A1      | A2
+               1950Q1 | 0.96284 | 0.5363
+               1950Q2 | 0.25145 | 0.31866
+               1950Q3 | 0.34447 | 0.4355
 
-                >> ts1 = dseries(rand(3,1),'1950Q2',{'A1'})
+               >> ts1 = dseries(rand(3,1),'1950Q2',{'A1'})
 
-                ts1 is a dseries object:
+               ts1 is a dseries object:
 
-                       | A1
-                1950Q2 | 0.70023
-                1950Q3 | 0.3958
-                1950Q4 | 0.084905
+                      | A1
+               1950Q2 | 0.40161
+               1950Q3 | 0.81763
+               1950Q4 | 0.97769
 
-                >> merge(ts0,ts1)
+               >> merge(ts0,ts1)
 
-                ans is a dseries object:
+               ans is a dseries object:
 
-                       | A1       | A2
-                1950Q1 | NaN      | 0.92477
-                1950Q2 | 0.70023  | 0.64208
-                1950Q3 | 0.3958   | 0.1045
-                1950Q4 | 0.084905 | NaN
+                      | A1      | A2
+               1950Q1 | 0.96284 | 0.5363
+               1950Q2 | 0.40161 | 0.31866
+               1950Q3 | 0.81763 | 0.4355
+               1950Q4 | 0.97769 | NaN
 
                 >> merge(ts1,ts0)
 
                 ans is a dseries object:
 
-                       | A1       | A2
-                1950Q1 | 0.42448  | 0.92477
-                1950Q2 | 0.60726  | 0.64208
-                1950Q3 | 0.070764 | 0.1045
-                1950Q4 | NaN      | NaN
-
-
-    .. dseriesmethod:: C = minus(A, B)
-
-        |br| Overloads the ``minus`` (``-``) operator for ``dseries``
-        objects, element by element subtraction. If both ``A`` and
-        ``B`` are ``dseries`` objects, they do not need to be defined
-        over the same time ranges. If ``A`` and ``B`` are ``dseries``
-        objects with :math:`T_A` and :math:`T_B` observations and
-        :math:`N_A` and :math:`N_B` variables, then :math:`N_A` must
-        be equal to :math:`N_B` or :math:`1` and :math:`N_B` must be
-        equal to :math:`N_A` or :math:`1`. If :math:`T_A=T_B`,
-        ``isequal(A.init,B.init)`` returns ``1`` and :math:`N_A=N_B`,
-        then the ``minus`` operator will compute for each couple
-        :math:`(t,n)`, with :math:`1\le t\le T_A` and :math:`1\le n\le
-        N_A`, ``C.data(t,n)=A.data(t,n)-B.data(t,n)``. If :math:`N_B`
-        is equal to :math:`1` and :math:`N_A>1`, the smaller
-        ``dseries`` object (``B``) is “broadcast” across the larger
-        ``dseries`` (``A``) so that they have compatible shapes, the
-        ``minus`` operator will subtract the variable defined in ``B``
-        from each variable in ``A``. If ``B`` is a double scalar, then
-        the method ``minus`` will subtract ``B`` from all the
+                      | A1      | A2
+               1950Q1 | 0.96284 | 0.5363
+               1950Q2 | 0.25145 | 0.31866
+               1950Q3 | 0.34447 | 0.4355
+               1950Q4 | 0.97769 | NaN
+
+
+    .. dseriesmethod:: C = minus (A, B)
+
+        |br| Overloads the MATLAB/Octave ``minus`` (``-``) operator
+        for ``dseries`` objects, element by element subtraction. If
+        both ``A`` and ``B`` are ``dseries`` objects, they do not need
+        to be defined over the same time ranges. If ``A`` and ``B``
+        are ``dseries`` objects with :math:`T_A` and :math:`T_B`
+        observations and :math:`N_A` and :math:`N_B` variables, then
+        :math:`N_A` must be equal to :math:`N_B` or :math:`1` and
+        :math:`N_B` must be equal to :math:`N_A` or :math:`1`. If
+        :math:`T_A=T_B`, ``isequal(A.init,B.init)`` returns ``1`` and
+        :math:`N_A=N_B`, then the ``minus`` operator will compute for
+        each couple :math:`(t,n)`, with :math:`1\le t\le T_A` and
+        :math:`1\le n\le N_A`,
+        ``C.data(t,n)=A.data(t,n)-B.data(t,n)``. If :math:`N_B` is
+        equal to :math:`1` and :math:`N_A>1`, the smaller ``dseries``
+        object (``B``) is “broadcast” across the larger ``dseries``
+        (``A``) so that they have compatible shapes, the ``minus``
+        operator will subtract the variable defined in ``B`` from each
+        variable in ``A``. If ``B`` is a double scalar, then the
+        method ``minus`` will subtract ``B`` from all the
         observations/variables in ``A``. If ``B`` is a row vector of
         length :math:`N_A`, then the ``minus`` method will subtract
         ``B(i)`` from all the observations of variable ``i``, for
@@ -1920,10 +2138,10 @@ The dseries class
 
                 ans is a dseries object:
 
-                   | minus(Variable_1,Variable_2) | minus(Variable_2,Variable_2)
-                1Y | -0.48853                     | 0
-                2Y | -0.50535                     | 0
-                3Y | -0.32063                     | 0
+                   | Variable_1 | Variable_2
+                1Y | -0.48853   | 0
+                2Y | -0.50535   | 0
+                3Y | -0.32063   | 0
 
                 >> ts1
 
@@ -1938,7 +2156,7 @@ The dseries class
 
                 ans is a dseries object:
 
-                   | minus(Variable_2,0.703)
+                   | Variable_2
                 1Y | 0
                 2Y | 0.051148
                 3Y | -0.15572
@@ -1947,15 +2165,15 @@ The dseries class
 
                 ans is a dseries object:
 
-                   | minus(0.703,Variable_2)
+                   | Variable_2
                 1Y | 0
                 2Y | -0.051148
                 3Y | 0.15572
 
 
-    .. dseriesmethod:: C = mpower(A, B)
+    .. dseriesmethod:: C = mpower (A, B)
 
-        |br| Overloads the ``mpower`` (``^``) operator for ``dseries``
+        |br| Overloads the MATLAB/Octave ``mpower`` (``^``) operator for ``dseries``
         objects and computes element-by-element power. ``A`` is a
         ``dseries`` object with ``N`` variables and ``T``
         observations. If ``B`` is a real scalar, then ``mpower(A,B)``
@@ -1974,7 +2192,7 @@ The dseries class
 
                 ts1 is a dseries object:
 
-                   | power(Variable_1,2)
+                   | Variable_1
                 1Y | 1
                 2Y | 4
                 3Y | 9
@@ -1983,15 +2201,15 @@ The dseries class
 
                 ts2 is a dseries object:
 
-                   | power(Variable_1,Variable_1)
+                   | Variable_1
                 1Y | 1
                 2Y | 4
                 3Y | 27
 
 
-    .. dseriesmethod:: C = mrdivide(A, B)
+    .. dseriesmethod:: C = mrdivide (A, B)
 
-        |br| Overloads the ``mrdivide`` (``/``) operator for
+        |br| Overloads the MATLAB/Octave ``mrdivide`` (``/``) operator for
         ``dseries`` objects, element by element division (like the
         ``./`` MATLAB/Octave operator). If both ``A`` and ``B`` are
         ``dseries`` objects, they do not need to be defined over the
@@ -2035,43 +2253,53 @@ The dseries class
 
                 ans is a dseries object:
 
-                   | divide(Variable_1,Variable_2) | divide(Variable_2,Variable_2)
-                1Y | 0.80745                       | 1
-                2Y | 4.2969                        | 1
-                3Y | 0.59235                       | 1
-
-
-    .. dseriesmethod:: C = mtimes(A, B)
-
-        |br| Overloads the ``mtimes`` (``*``) operator for ``dseries``
-        objects and the Hadammard product (the .* MATLAB/Octave
-        operator). If both ``A`` and ``B`` are ``dseries`` objects,
-        they do not need to be defined over the same time ranges. If
-        ``A`` and ``B`` are ``dseries`` objects with :math:`T_A` and
-        :math:`_B` observations and :math:`N_A` and :math:`N_B`
-        variables, then :math:`N_A` must be equal to :math:`N_B` or
-        :math:`1` and :math:`N_B` must be equal to :math:`N_A` or
-        :math:`1`. If :math:`T_A=T_B`, ``isequal(A.init,B.init)``
-        returns ``1`` and :math:`N_A=N_B`, then the ``mtimes``
-        operator will compute for each couple :math:`(t,n)`, with
-        :math:`1\le t\le T_A` and :math:`1\le n\le N_A`,
-        ``C.data(t,n)=A.data(t,n)*B.data(t,n)``. If :math:`N_B` is
-        equal to :math:`1` and :math:`N_A>1`, the smaller ``dseries``
-        object (``B``) is “broadcast” across the larger ``dseries``
-        (``A``) so that they have compatible shapes, ``mtimes``
-        operator will multiply each variable defined in ``A`` by the
-        variable in ``B``, observation per observation. If ``B`` is a
-        double scalar, then the method ``mtimes`` will multiply all
-        the observations/variables in ``A`` by ``B``. If ``B`` is a
-        row vector of length :math:`N_A`, then the ``mtimes`` method
-        will multiply all the observations of variable ``i`` by
-        ``B(i)``, for :math:`i=1,...,N_A`. If ``B`` is a column vector
-        of length :math:`T_A`, then the ``mtimes`` method will perform
-        a multiplication of all the variables by ``B``, element by
+                   | Variable_1 | Variable_2
+                1Y | 0.80745    | 1
+                2Y | 4.2969     | 1
+                3Y | 0.59235    | 1
+
+
+    .. dseriesmethod:: C = mtimes (A, B)
+
+        |br| Overloads the MATLAB/Octave ``mtimes`` (``*``) operator
+        for ``dseries`` objects and the Hadammard product (the .*
+        MATLAB/Octave operator). If both ``A`` and ``B`` are
+        ``dseries`` objects, they do not need to be defined over the
+        same time ranges. If ``A`` and ``B`` are ``dseries`` objects
+        with :math:`T_A` and :math:`_B` observations and :math:`N_A`
+        and :math:`N_B` variables, then :math:`N_A` must be equal to
+        :math:`N_B` or :math:`1` and :math:`N_B` must be equal to
+        :math:`N_A` or :math:`1`. If :math:`T_A=T_B`,
+        ``isequal(A.init,B.init)`` returns ``1`` and :math:`N_A=N_B`,
+        then the ``mtimes`` operator will compute for each couple
+        :math:`(t,n)`, with :math:`1\le t\le T_A` and :math:`1\le n\le
+        N_A`, ``C.data(t,n)=A.data(t,n)*B.data(t,n)``. If :math:`N_B`
+        is equal to :math:`1` and :math:`N_A>1`, the smaller
+        ``dseries`` object (``B``) is “broadcast” across the larger
+        ``dseries`` (``A``) so that they have compatible shapes,
+        ``mtimes`` operator will multiply each variable defined in
+        ``A`` by the variable in ``B``, observation per
+        observation. If ``B`` is a double scalar, then the method
+        ``mtimes`` will multiply all the observations/variables in
+        ``A`` by ``B``. If ``B`` is a row vector of length
+        :math:`N_A`, then the ``mtimes`` method will multiply all the
+        observations of variable ``i`` by ``B(i)``, for
+        :math:`i=1,...,N_A`. If ``B`` is a column vector of length
+        :math:`T_A`, then the ``mtimes`` method will perform a
+        multiplication of all the variables by ``B``, element by
         element.
 
 
-    .. dseriesmethod:: C = ne(A, B)
+    .. dseriesmethod:: B = nanmean (A[, geometric])
+
+        |br| Overloads the MATLAB/Octave ``nanmean`` function for
+        ``dseries`` objects. Returns the mean of each variable in
+        ``dseries`` object ``A`` ignoring the NaN values. If the
+        second argument is ``true`` the geometric mean is computed,
+        otherwise (default) the arithmetic mean is reported.
+
+
+    .. dseriesmethod:: C = ne (A, B)
 
         |br| Overloads the MATLAB/Octave ``ne`` (not equal, ``~=``)
         operator. ``dseries`` objects ``A`` and ``B`` must have the
@@ -2091,12 +2319,14 @@ The dseries class
 
                 ans =
 
-                     0
-                     1
-                     0
+                  3x1 logical array
+
+                   0
+                   1
+                   0
 
 
-    .. dseriesmethod:: B = nobs(A)
+    .. dseriesmethod:: B = nobs (A)
 
         |br| Returns the number of observations in ``dseries`` object
         ``A``.
@@ -2113,10 +2343,32 @@ The dseries class
                     10
 
 
-    .. dseriesmethod:: h = plot(A)
-                       h = plot(A, B)
-                       h = plot(A[, ...])
-                       h = plot(A, B[, ...])
+    .. dseriesmethod:: B = onesidedhpcycle (A[, lambda[, init]])
+                       onesidedhpcycle_ (A[, lambda[, init]])
+
+        |br| Extracts the cycle component from a ``dseries`` ``A``
+        object using a one sided HP filter (with a Kalman filter) and
+        returns a ``dseries`` object, ``B``. The default value for
+        ``lambda``, the smoothing parameter, is ``1600``. By default,
+        if ``ìnit`` is not provided, the initial value is based on the
+        first two observations.
+
+
+    .. dseriesmethod:: B = onesidedhptrend (A[, lambda[, init]])
+                       onesidedhptrend_ (A[, lambda[, init]])
+
+        |br| Extracts the trend component from a ``dseries`` ``A``
+        object using a one sided HP filter (with a Kalman filter) and
+        returns a ``dseries`` object, ``B``. The default value for
+        ``lambda``, the smoothing parameter, is ``1600``. By default,
+        if ``ìnit`` is not provided, the initial value is based on the
+        first two observations.
+
+
+    .. dseriesmethod:: h = plot (A)
+                       h = plot (A, B)
+                       h = plot (A[, ...])
+                       h = plot (A, B[, ...])
 
         |br| Overloads MATLAB/Octave’s ``plot`` function for
         ``dseries`` objects. Returns a MATLAB/Octave plot handle, that
@@ -2176,35 +2428,36 @@ The dseries class
                 >> set(h(2),'+r');
 
 
-    .. dseriesmethod:: C = plus(A, B)
-
-        |br| Overloads the ``plus`` (``+``) operator for ``dseries``
-        objects, element by element addition. If both ``A`` and ``B``
-        are ``dseries`` objects, they do not need to be defined over
-        the same time ranges. If ``A`` and ``B`` are ``dseries``
-        objects with :math:`T_A` and :math:`T_B` observations and
-        :math:`N_A` and :math:`N_B` variables, then :math:`N_A` must
-        be equal to :math:`N_B` or :math:`1` and :math:`N_B` must be
-        equal to :math:`N_A` or :math:`1`. If :math:`T_A=T_B`,
-        ``isequal(A.init,B.init)`` returns ``1`` and :math:`N_A=N_B`,
-        then the ``plus`` operator will compute for each couple
-        :math:`(t,n)`, with :math:`1\le t\le T_A` and :math:`1\le n\le
-        N_A`, ``C.data(t,n)=A.data(t,n)+B.data(t,n)``. If :math:`N_B`
-        is equal to :math:`1` and :math:`N_A>1`, the smaller
-        ``dseries`` object (``B``) is “broadcast” across the larger
-        ``dseries`` (``A``) so that they have compatible shapes, the
-        plus operator will add the variable defined in ``B`` to each
-        variable in ``A``. If ``B`` is a double scalar, then the
-        method ``plus`` will add ``B`` to all the
-        observations/variables in ``A``. If ``B`` is a row vector of
-        length :math:`N_A`, then the ``plus`` method will add ``B(i)``
-        to all the observations of variable ``i``, for
-        :math:`i=1,...,N_A`. If ``B`` is a column vector of length
-        :math:`T_A`, then the ``plus`` method will add ``B`` to all
-        the variables.
+    .. dseriesmethod:: C = plus (A, B)
+
+        |br| Overloads the MATLAB/Octave ``plus`` (``+``) operator for
+        ``dseries`` objects, element by element addition. If both
+        ``A`` and ``B`` are ``dseries`` objects, they do not need to
+        be defined over the same time ranges. If ``A`` and ``B`` are
+        ``dseries`` objects with :math:`T_A` and :math:`T_B`
+        observations and :math:`N_A` and :math:`N_B` variables, then
+        :math:`N_A` must be equal to :math:`N_B` or :math:`1` and
+        :math:`N_B` must be equal to :math:`N_A` or :math:`1`. If
+        :math:`T_A=T_B`, ``isequal(A.init,B.init)`` returns ``1`` and
+        :math:`N_A=N_B`, then the ``plus`` operator will compute for
+        each couple :math:`(t,n)`, with :math:`1\le t\le T_A` and
+        :math:`1\le n\le N_A`,
+        ``C.data(t,n)=A.data(t,n)+B.data(t,n)``. If :math:`N_B` is
+        equal to :math:`1` and :math:`N_A>1`, the smaller ``dseries``
+        object (``B``) is “broadcast” across the larger ``dseries``
+        (``A``) so that they have compatible shapes, the plus operator
+        will add the variable defined in ``B`` to each variable in
+        ``A``. If ``B`` is a double scalar, then the method ``plus``
+        will add ``B`` to all the observations/variables in ``A``. If
+        ``B`` is a row vector of length :math:`N_A`, then the ``plus``
+        method will add ``B(i)`` to all the observations of variable
+        ``i``, for :math:`i=1,...,N_A`. If ``B`` is a column vector of
+        length :math:`T_A`, then the ``plus`` method will add ``B`` to
+        all the variables.
 
 
-    .. dseriesmethod:: C = pop(A[, B])
+    .. dseriesmethod:: C = pop (A[, B])
+                       pop_ (A[, B])
 
         |br| Removes variable ``B`` from ``dseries`` object ``A``. By
         default, if the second argument is not provided, the last
@@ -2225,8 +2478,10 @@ The dseries class
                 3Y | 1          | 1
 
 
-    .. dseriesmethod:: B = qdiff(A)
-                       B = qgrowth(A)
+    .. dseriesmethod:: B = qdiff (A)
+                       B = qgrowth (A)
+                       qdiff_ (A)
+                       qgrowth_ (A)
 
         |br| Computes quarterly differences or growth rates.
 
@@ -2239,7 +2494,7 @@ The dseries class
 
                 ts1 is a dseries object:
 
-                       | qdiff(Variable_1)
+                       | Variable_1
                 1950Q1 | NaN
                 1950Q2 | 1
                 1950Q3 | 1
@@ -2250,7 +2505,7 @@ The dseries class
 
                 ts1 is a dseries object:
 
-                        | qdiff(Variable_1)
+                        | Variable_1
                 1950M1  | NaN
                 1950M2  | NaN
                 1950M3  | NaN
@@ -2259,7 +2514,8 @@ The dseries class
                 1950M6  | 3
 
 
-    .. dseriesmethod:: C = remove(A, B)
+    .. dseriesmethod:: C = remove (A, B)
+                       remove_ (A, B)
 
         |br| Alias for the ``pop`` method with two arguments. Removes
         variable ``B`` from ``dseries`` object ``A``.
@@ -2292,10 +2548,13 @@ The dseries class
             implicit loops can.
 
 
-    .. dseriesmethod:: B = rename(A,oldname,newname)
+    .. dseriesmethod:: B = rename (A, oldname, newname)
+                       rename_ (A, oldname, newname)
 
         |br| Rename variable ``oldname`` to ``newname`` in ``dseries``
-        object ``A``. Returns a ``dseries`` object.``
+        object ``A``. Returns a ``dseries`` object. If more than one
+        variable needs to be renamed, it is possible to pass cells of
+        char arrays as second and third arguments.
 
         *Example*
 
@@ -2311,12 +2570,13 @@ The dseries class
                 2Y | 1       | 1
 
 
-    .. dseriesmethod:: C = rename(A,newname)
+    .. dseriesmethod:: C = rename (A, newname)
+                       rename_ (A, newname)
 
         |br| Replace the names in ``A`` with those passed in the cell
         string array ``newname``. ``newname`` must have the same
-        number of cells as ``A`` has ``dseries``. Returns a
-        ``dseries`` object.
+        number of elements as ``dseries`` object ``A`` has
+        variables. Returns a ``dseries`` object.
 
         *Example*
 
@@ -2332,12 +2592,12 @@ The dseries class
                 2Y | 1          | 1     | 1
 
 
-    .. dseriesmethod:: save(A, basename[, format])
+    .. dseriesmethod:: save (A, basename[, format])
 
         |br| Overloads the MATLAB/Octave ``save`` function and saves
-        ``dseries`` object ``A`` to disk. Possible formats are ``csv``
+        ``dseries`` object ``A`` to disk. Possible formats are ``mat``
         (this is the default), ``m`` (MATLAB/Octave script), and
-        ``mat`` (MATLAB binary data file). The name of the file
+        ``csv`` (MATLAB binary data file). The name of the file
         without extension is specified by ``basename``.
 
         *Example*
@@ -2345,7 +2605,7 @@ The dseries class
             ::
 
                 >> ts0 = dseries(ones(2,2));
-                >> ts0.save('ts0');
+                >> ts0.save('ts0', 'csv');
 
             The last command will create a file ts0.csv with the
             following content::
@@ -2367,6 +2627,8 @@ The dseries class
 
                 NAMES__ = {'Variable_1'; 'Variable_2'};
                 TEX__ = {'Variable_{1}'; 'Variable_{2}'};
+                OPS__ = {};
+                TAGS__ = struct();
 
                 Variable_1 = [
                               1
@@ -2426,8 +2688,34 @@ The dseries class
                      1     3
 
 
-    .. dseriesmethod:: B = tex_rename(A, name, newtexname)
-                       B = tex_rename(A, newtexname)
+    .. dseriesmethod:: B = std (A[, geometric])
+
+        |br| Overloads the MATLAB/Octave ``std`` function for
+        ``dseries`` objects. Returns the standard deviation of each
+        variable in ``dseries`` object ``A``. If the second argument
+        is ``true`` the geometric standard deviation is computed
+        (default value of the second argument is ``false``).
+
+
+    .. dseriesmethod:: A = tag (A, a[, b, c])
+
+        |br| Add a tag to a variable in ``dseries`` object ``A``.
+
+        *Example*
+
+            ::
+
+               >> ts = dseries(randn(10, 3));
+               >> tag(ts, 'type');             % Define a tag name.
+               >> tag(ts, 'type', 'Variable_1', 'Stock');
+               >> tag(ts, 'type', 'Variable_2', 'Flow');
+               >> tag(ts, 'type', 'Variable_3', 'Stock');
+
+
+    .. dseriesmethod:: B = tex_rename (A, name, newtexname)
+                       B = tex_rename (A, newtexname)
+                       tex_rename_ (A, name, newtexname)
+                       tex_rename_ (A, newtexname)
 
         |br| Redefines the tex name of variable ``name`` to
         ``newtexname`` in ``dseries`` object ``A``. Returns a
@@ -2459,7 +2747,7 @@ The dseries class
 
                 ts1 is a dseries object:
 
-                   | -Variable_1
+                   | Variable_1
                 1Y | -1
 
 
@@ -2490,7 +2778,7 @@ The dseries class
                 1950Q4 | 0.11171  | 0.67865
 
 
-    .. dseriesmethod:: B = vobs(A)
+    .. dseriesmethod:: B = vobs (A)
 
         |br| Returns the number of variables in ``dseries`` object
         ``A``.
@@ -2507,7 +2795,9 @@ The dseries class
                     2
 
 
-.. dseriesmethod:: B = ydiff(A)
-                   B = ygrowth(A)
+.. dseriesmethod:: B = ydiff (A)
+                   B = ygrowth (A)
+                   ydiff_ (A)
+                   ygrowth_ (A)
 
         |br| Computes yearly differences or growth rates.