From 5aab3b5aa95beca09a09d8aae57a51e4e977c1df Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Adjemian=20=28Charybdis=29?=
 <stepan@adjemian.eu>
Date: Thu, 8 Oct 2020 23:32:28 +0200
Subject: [PATCH] Added routines to convert dseries to lower frequencies.

 - daily to monthly
 - daily or monthly to quaterly

Conversion to annual frequency is not yet implemented, but this is
enough to close #1.

Four time aggregation methods are available:

 - 'sum' (aggregation of flow variables),
 - 'end-of-period' (aggregation of stock variables),
 - 'arithmetic-average', (aggregation of rates), and
 - 'geometric-average' (aggregation of ratios)

Examples:

>> a = dseries(rand(2*365+1,1), '2000-01-01')

a is a dseries object:

           | Variable_1
2000-01-01 | 0.26032
2000-01-02 | 0.56927
2000-01-03 | 0.24877
2000-01-04 | 0.3193
2000-01-05 | 0.9108
2000-01-06 | 0.88522
2000-01-07 | 0.79459
2000-01-08 | 0.92581
2000-01-09 | 0.17884
2000-01-10 | 0.51754
           |
2001-12-21 | 0.035016
2001-12-22 | 0.03975
2001-12-23 | 0.98857
2001-12-24 | 0.68618
2001-12-25 | 0.37669
2001-12-26 | 0.50432
2001-12-27 | 0.7635
2001-12-28 | 0.048875
2001-12-29 | 0.72593
2001-12-30 | 0.70133
2001-12-31 | 0.45889

>> dseries2Q(a, 'sum')

ans is a dseries object:

       | Variable_1
2000Q1 | 46.6096
2000Q2 | 45.1456
2000Q3 | 45.0172
2000Q4 | 47.6164
2001Q1 | 49.9093
2001Q2 | 46.2334
2001Q3 | 47.4685
2001Q4 | 47.1611

>> dseries2M(a, 'sum')

ans is a dseries object:

        | Variable_1
2000M1  | 16.7499
2000M2  | 15.4311
2000M3  | 14.4286
2000M4  | 14.3583
2000M5  | 15.7072
2000M6  | 15.0801
2000M7  | 14.2406
2000M8  | 15.9174
2000M9  | 14.8592
2000M10 | 17.5525
2000M11 | 13.3236
2000M12 | 16.7404
2001M1  | 18.9293
2001M2  | 16.1282
2001M3  | 14.8518
2001M4  | 16.9843
2001M5  | 14.5768
2001M6  | 14.6723
2001M7  | 14.9768
2001M8  | 15.3968
2001M9  | 17.0949
2001M10 | 15.0376
2001M11 | 16.6848
2001M12 | 15.4388

>> dseries2Q(dseries2M(a, 'sum'), 'sum')

ans is a dseries object:

       | Variable_1
2000Q1 | 46.6096
2000Q2 | 45.1456
2000Q3 | 45.0172
2000Q4 | 47.6164
2001Q1 | 49.9093
2001Q2 | 46.2334
2001Q3 | 47.4685
2001Q4 | 47.1611
---
 src/utilities/convert/dseries2M.m | 102 ++++++++++++++++++++++++
 src/utilities/convert/dseries2Q.m | 127 ++++++++++++++++++++++++++++++
 2 files changed, 229 insertions(+)
 create mode 100644 src/utilities/convert/dseries2M.m
 create mode 100644 src/utilities/convert/dseries2Q.m

diff --git a/src/utilities/convert/dseries2M.m b/src/utilities/convert/dseries2M.m
new file mode 100644
index 0000000..9b3768f
--- /dev/null
+++ b/src/utilities/convert/dseries2M.m
@@ -0,0 +1,102 @@
+function ts = dseries2M(ds, method)
+
+% INPUTS
+% - ds       [dseries]    daily frequency object with n elements.
+% - method   [char]       1×m array, method used for the time aggregation, possible values are:
+%                            - 'arithmetic-average' (for rates),
+%                            - 'geometric-average' (for factors),
+%                            - 'sum' (for flow variables), and
+%                            - 'end-of-period' (for stock variables).
+%
+% OUTPUTS
+% - ts       [dseries]    monthly frequency object with m<n elements.
+
+% Copyright © 2020 Dynare Team
+%
+% This code 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 dates submodule 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 <http://www.gnu.org/licenses/>.
+
+    dm = dates2M(ds.dates);
+    dM = unique(dm); % dM.ndat will be the maximum number of periods in ts
+
+    tdata = NaN(dM.ndat, ds.vobs);
+
+    switch method
+      case 'arithmetic-average'
+        for i=1:dM.ndat
+            idm = find(dm==dM(i));
+            tdata(i,:) = mean(ds.data(idm,:));
+            if ~isalldays(dM(i), idm)
+                warning('Not all days are available in month %s.', date2string(dM(i)))
+            end
+        end
+      case 'geometric-average'
+        for i=1:dM.ndat
+            idm = find(dm==dM(i));
+            tdata(i,:) = prod(ds.data(idm,:), 1).^(1/length(idm));
+            if ~isalldays(dM(i), idm)
+                warning('Not all days are available in month %s.', date2string(dM(i)))
+            end
+        end
+      case 'sum'
+        for i=1:dM.ndat
+            idm = find(dm==dM(i));
+            tdata(i,:) = sum(ds.data(idm,:), 1);
+            if ~isalldays(dM(i), idm)
+                warning('Not all days are available in month %s.', date2string(dM(i)))
+            end
+        end
+      case 'end-of-period'
+        for i=1:dM.ndat
+            idm = find(dm==dM(i));
+            lastday = datevec(ds.dates.time(idm(end)));
+            tdata(i,:) = ds.data(idm(end),:);
+            if lastday(3)<expectedlastday(dM(i))
+                warning('The last available day is not the last observed day of %s.', date2string(dM(i)))
+            end
+        end
+      otherwise
+        error('Unknow method.')
+    end
+
+    ts = dseries(tdata, dM, ds.name, ds.tex);
+
+end
+
+function eld = expectedlastday(m)
+    if ismember(m.time(2), [4,6,9,11])
+        eld = 30;
+    elseif ismember(m.time(2), [1,3,5,7,8,10,12])
+        eld = 31;
+    else % February
+        if isleapyear(m.time(1))
+            eld = 29;
+        else
+            eld = 28;
+        end
+    end
+end
+
+function iad = isalldays(m, id)
+    if ismember(m.time(2), [4,6,9,11])
+        iad = length(id)==30;
+    elseif ismember(m.time(2), [1,3,5,7,8,10,12])
+        iad = length(id)==31;
+    else % February
+        if isleapyear(m.time(1))
+            iad = length(id)==29;
+        else
+            iad = length(id)==28;
+        end
+    end
+end
\ No newline at end of file
diff --git a/src/utilities/convert/dseries2Q.m b/src/utilities/convert/dseries2Q.m
new file mode 100644
index 0000000..0e05bac
--- /dev/null
+++ b/src/utilities/convert/dseries2Q.m
@@ -0,0 +1,127 @@
+function ts = dseries2Q(ds, method)
+
+% INPUTS
+% - ds       [dseries]    daily or monthly frequency object with n elements.
+% - method   [char]       1×m array, method used for the time aggregation, possible values are:
+%                            - 'arithmetic-average' (for rates),
+%                            - 'geometric-average' (for factors),
+%                            - 'sum' (for flow variables), and
+%                            - 'end-of-period' (for stock variables).
+%
+% OUTPUTS
+% - ts       [dseries]    quaterly frequency object with m<n elements.
+
+% Copyright © 2020 Dynare Team
+%
+% This code 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 dates submodule 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 <http://www.gnu.org/licenses/>.
+
+    dq = dates2Q(ds.dates);
+    dQ = unique(dq); % dQ.ndat will be the maximum number of periods in ts
+
+    tdata = NaN(dQ.ndat, ds.vobs);
+
+    switch method
+      case 'arithmetic-average'
+        for i=1:dQ.ndat
+            idq = find(dq==dQ(i));
+            tdata(i,:) = mean(ds.data(idq,:));
+            if ds.dates.freq==365 && ~isalldays(dQ(i), idq)
+                warning('Not all days are available in %s.', date2string(dQ(i)))
+            end
+            if ds.dates.freq==12 && ~isequal(length(idq), 3)
+                warning('Not all months are available in %s.', date2string(dQ(i)))
+            end
+        end
+      case 'geometric-average'
+        for i=1:dQ.ndat
+            idq = find(dq==dQ(i));
+            tdata(i,:) = prod(ds.data(idq,:), 1).^(1/length(idq));
+            if ds.dates.freq==365 && ~isalldays(dQ(i), idq)
+                warning('Not all days are available in %s.', date2string(dQ(i)))
+            end
+            if ds.dates.freq==12 && ~isequal(length(idq), 3)
+                warning('Not all months are available in %s.', date2string(dQ(i)))
+            end
+        end
+      case 'sum'
+        for i=1:dQ.ndat
+            idq = find(dq==dQ(i));
+            tdata(i,:) = sum(ds.data(idq,:), 1);
+            if ds.dates.freq==365 && ~isalldays(dQ(i), idq)
+                warning('Not all days are available in %s.', date2string(dQ(i)))
+            end
+            if ds.dates.freq==12 && ~isequal(length(idq), 3)
+                warning('Not all months are available in %s.', date2string(dQ(i)))
+            end
+        end
+      case 'end-of-period'
+        for i=1:dQ.ndat
+            idq = find(dq==dQ(i));
+            tdata(i,:) = ds.data(idq(end),:);
+            if ds.dates.freq==365
+                lastday = datevec(ds.dates.time(idq(end)));
+                if lastday(3)<expectedlastday(dQ(i))
+                    warning('The last available day is not the last observed day of %s.', date2string(dQ(i)))
+                end
+            end
+            if ds.dates.freq==12
+                lastmonth = ds.dates.time(idq(end), 2);
+                if lastmonth<expectedlastmonth(dQ(i))
+                    warning('The last available day is not the last observed day of %s.', date2string(dQ(i)))
+                end
+            end
+        end
+      otherwise
+        error('Unknow method.')
+    end
+
+    ts = dseries(tdata, dQ, ds.name, ds.tex);
+
+end
+
+function eld = expectedlastday(q)
+    if ismember(q.time(2), [1,4])
+        eld = 31; % last day of March or December
+    else
+        eld = 30; % last day of June or September
+    end
+end
+
+function elm = expectedlastmonth(q)
+    switch q.time(2)
+      case 1
+        elm = 3;
+      case 2
+        elm = 6;
+      case 3
+        elm = 9;
+      case 4
+        elm = 12;
+      otherwise
+    end
+end
+
+function iad = isalldays(m, id)
+    if ismember(m.time(2), [3, 4])
+        iad = length(id)==92;
+    elseif m.time(2)==2
+        iad = length(id)==91;
+    else % contains February
+        if isleapyear(m.time(1))
+            iad = length(id)==91;
+        else
+            iad = length(id)==90;
+        end
+    end
+end
\ No newline at end of file
-- 
GitLab