From a3c6f14776f299b7eeee22c4dec2c63352e68302 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Adjemian=20=28Charybdis=29?=
 <stephane.adjemian@univ-lemans.fr>
Date: Fri, 21 Oct 2016 15:30:14 +0200
Subject: [PATCH] Changed behaviour of cumsum and cumprod...

When NaNs are present in the dataset. Do not return an array filled with
NaNs if only the first observations are NaNs.
---
 src/@dseries/cumprod.m              | 50 +++++++++++++++++++++++----
 src/@dseries/cumsum.m               | 52 +++++++++++++++++++++++++----
 src/initialize_dseries_toolbox.m    |  1 +
 src/utilities/cumulate/cumprodnan.m | 24 +++++++++++++
 src/utilities/cumulate/cumsumnan.m  | 24 +++++++++++++
 5 files changed, 137 insertions(+), 14 deletions(-)
 create mode 100644 src/utilities/cumulate/cumprodnan.m
 create mode 100644 src/utilities/cumulate/cumsumnan.m

diff --git a/src/@dseries/cumprod.m b/src/@dseries/cumprod.m
index 0cd0247..4e79bae 100644
--- a/src/@dseries/cumprod.m
+++ b/src/@dseries/cumprod.m
@@ -27,19 +27,40 @@ function B = cumprod(varargin) % --*-- Unitary tests --*--
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
 
-% Get indices of the columns without NaNs
-idx = find(~any(isnan(varargin{1}.data)));
+% Get the firt observation number where all the variables are observed (ie without NaNs)
+idx = find(all(~isnan(varargin{1}.data), 2),1);
+if isempty(idx)
+    idx = 0;
+end
 
-if ~isequal(idx(:),transpose(1:vobs(varargin{1})))
-    warning('dseries::cumprod: The cumulated product is not computed for some variables because they have NaNs!')
+% Is the first period where variables are observed common?
+common_first_period_witout_nan = true;
+if ~idx
+    if any(~isnan(varargin{1}.data(:)))
+        common_first_period_witout_nan = false;
+    end
+else
+    if idx>1
+        if any(any(~isnan(varargin{1}.data(1:idx-1,:))))
+            common_first_period_witout_nan = false;
+        end
+    end
 end
 
 switch nargin
-    case 1
+  case 1
       % Initialize the output.
       B = varargin{1};
       % Perform the cumulated sum
-      B.data(:,idx) = cumprod(B.data(:,idx));
+      if isequal(idx, 1)
+          B.data = cumprod(B.data);
+      else
+          if common_first_period_witout_nan
+              B.data(idx:end,:) = cumprod(B.data(idx:end,:));
+          else
+              B.data = cumprodnan(B.data);
+          end
+      end
       % Change the name of the variables
       for i=1:vobs(B)
           B.name(i) = {['cumprod(' B.name{i} ')']};
@@ -221,4 +242,19 @@ end
 %$ % Check the results.
 %$ t(1) = dassert(ts3.data,ts4.data);
 %$ T = all(t);
-%@eof:5
\ No newline at end of file
+%@eof:5
+
+%@test:6
+%$ % Define a data set.
+%$ A = [NaN, NaN; 1 NaN; 1 1; 1 1; 1 NaN];
+%$
+%$ % Instantiate a time series object.
+%$ ts0 = dseries(A);
+%$
+%$ % Call the tested method.
+%$ ts0.cumprod();
+%$
+%$ % Check the results.
+%$ t(1) = dassert(ts0.data, A);
+%$ T = all(t);
+%@eof:6
\ No newline at end of file
diff --git a/src/@dseries/cumsum.m b/src/@dseries/cumsum.m
index 7acaea3..6b6ed11 100644
--- a/src/@dseries/cumsum.m
+++ b/src/@dseries/cumsum.m
@@ -27,12 +27,24 @@ function B = cumsum(varargin) % --*-- Unitary tests --*--
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
 
-% Get indices of the columns without NaNs
-idx = find(~any(isnan(varargin{1}.data)));
-
+% Get the firt observation number where all the variables are observed (ie without NaNs)
+idx = find(all(~isnan(varargin{1}.data), 2),1);
+if isempty(idx)
+    idx = 0;
+end
 
-if ~isequal(idx(:),transpose(1:vobs(varargin{1})))
-    warning('dseries::cumsum: The cumulated sum is not computed for some variables because they have NaNs!')
+% Is the first period where variables are observed common?
+common_first_period_witout_nan = true;
+if ~idx
+    if any(~isnan(varargin{1}.data(:)))
+        common_first_period_witout_nan = false;
+    end
+else
+    if idx>1
+        if any(any(~isnan(varargin{1}.data(1:idx-1,:))))
+            common_first_period_witout_nan = false;
+        end
+    end
 end
 
 switch nargin
@@ -40,7 +52,15 @@ switch nargin
       % Initialize the output.
       B = varargin{1};
       % Perform the cumulated sum
-      B.data(:,idx) = cumsum(B.data(:,idx));
+      if isequal(idx, 1)
+          B.data = cumsum(B.data);
+      else
+          if common_first_period_witout_nan
+              B.data(idx:end,:) = cumsum(B.data(idx:end,:));
+          else
+              B.data = cumsumnan(B.data);
+          end
+      end
       % Change the name of the variables
       for i=1:vobs(B)
           B.name(i) = {['cumsum(' B.name{i} ')']};
@@ -209,4 +229,22 @@ end
 %$ % Check the results.
 %$ t(1) = dassert(ts3.data,ts4.data);
 %$ T = all(t);
-%@eof:5
\ No newline at end of file
+%@eof:5
+
+%@test:6
+%$ % Define a data set.
+%$ A = [NaN, NaN; 1 NaN; 1 1; 1 1; 1 NaN];
+%$
+%$ % Instantiate a time series object.
+%$ ts0 = dseries(A);
+%$
+%$ % Call the tested method.
+%$ ts0.cumsum();
+%$
+%$ % Expected results.
+%$ A = [NaN   NaN; 1   NaN; 2     1; 3     2; 4   NaN];
+%$
+%$ % Check the results.
+%$ t(1) = dassert(ts0.data, A);
+%$ T = all(t);
+%@eof:6
\ No newline at end of file
diff --git a/src/initialize_dseries_toolbox.m b/src/initialize_dseries_toolbox.m
index 60719e3..0d8ff39 100644
--- a/src/initialize_dseries_toolbox.m
+++ b/src/initialize_dseries_toolbox.m
@@ -23,6 +23,7 @@ addpath([dseries_src_root '/utilities/insert'])
 addpath([dseries_src_root '/utilities/file'])
 addpath([dseries_src_root '/utilities/from'])
 addpath([dseries_src_root '/utilities/variables'])
+addpath([dseries_src_root '/utilities/cumulate'])
 
 % Add missing routines if dynare is not in the path
 if ~exist('demean','file')
diff --git a/src/utilities/cumulate/cumprodnan.m b/src/utilities/cumulate/cumprodnan.m
new file mode 100644
index 0000000..f0d3a78
--- /dev/null
+++ b/src/utilities/cumulate/cumprodnan.m
@@ -0,0 +1,24 @@
+function A = cumprodnan(a)
+
+% Copyright (C) 2016 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 dseries 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/>.
+
+b = ~isnan(a);
+A = a;
+
+for i=1:size(a, 2)
+    c = find(b(:,i), 1);
+    A(c:end,i) = cumprod(a(c:end,i));
+end
\ No newline at end of file
diff --git a/src/utilities/cumulate/cumsumnan.m b/src/utilities/cumulate/cumsumnan.m
new file mode 100644
index 0000000..2ee14ae
--- /dev/null
+++ b/src/utilities/cumulate/cumsumnan.m
@@ -0,0 +1,24 @@
+function A = cumsumnan(a)
+
+% Copyright (C) 2016 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 dseries 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/>.
+
+b = ~isnan(a);
+A = a;
+
+for i=1:size(a, 2)
+    c = find(b(:,i), 1);
+    A(c:end,i) = cumsum(a(c:end,i));
+end
\ No newline at end of file
-- 
GitLab