← Back to branch summary

~ubuntu-branches/debian/sid/octave-tisean/sid

« back to all changes in this revision

Viewing changes to inst/timerev.m

  • Committer: Package Import Robot
  • Author(s): Rafael Laboissiere
  • Date: 2017-08-14 12:53:47 UTC
  • Revision ID: package-import@ubuntu.com-20170814125347-ju5owr4dggr53a2n
Tags: upstream-0.2.3
Import upstream version 0.2.3

Show diffs side-by-side

added added

removed removed

Lines of Context:
inst/timerev.m
 
1
## Copyright (C) 2015 Piotr Held
 
2
## Copyright (C) 2015 Juan Pablo Carbajal
 
3
##
 
4
## This file is part of Octave.
 
5
##
 
6
## Octave is free software; you can redistribute it and/or
 
7
## modify it under the terms of the GNU General Public
 
8
## License as published by the Free Software Foundation;
 
9
## either version 3 of the License, or (at your option) any
 
10
## later version.
 
11
##
 
12
## Octave is distributed in the hope that it will be useful,
 
13
## but WITHOUT ANY WARRANTY; without even the implied
 
14
## warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
 
15
## PURPOSE.  See the GNU General Public License for more
 
16
## details.
 
17
##
 
18
## You should have received a copy of the GNU General Public
 
19
## License along with Octave; see the file COPYING.  If not,
 
20
## see <http://www.gnu.org/licenses/>.
 
21
 
 
22
## -*- texinfo -*-
 
23
## @deftypefn{Function File} {@var{output} =} timerev (@var{S})
 
24
## @deftypefnx{Function File} {@var{output} =} timerev (@var{S}, @var{delay})
 
25
##
 
26
## Calculates time reversal assymetry statistic.
 
27
##
 
28
## Accomplishes this using the following equation applied to each component
 
29
## separately:
 
30
##
 
31
## @iftex
 
32
## @tex
 
33
## $$\frac{\sum (y_n-y_{n-d})^3}{\sum (y_n-y_{n-d})^2}$$
 
34
## @end tex
 
35
## @end iftex
 
36
## @ifnottex
 
37
## @example
 
38
##                 3
 
39
##  sum (y  - y   ) 
 
40
##        n    n-d
 
41
## ------------------
 
42
##                 2
 
43
##  sum (y  - y   )  
 
44
##        n    n-d
 
45
## @end example
 
46
## @end ifnottex
 
47
##
 
48
## @strong{Input}
 
49
##
 
50
## @table @var
 
51
## @item S
 
52
## This function always assumes that each time series is along the longer 
 
53
## dimension of matrix @var{S}. It also assumes that every dimension 
 
54
## (counting along the shorter dimension) of @var{S} is considered a 
 
55
## component of the time series.
 
56
## @item delay
 
57
## The delay for the statistic ('d' in the equation above) [default = 1].
 
58
## @end table
 
59
##
 
60
## @strong{Output}
 
61
##
 
62
## The output is the calculated time reversal asymmetry statistic. It is calculated
 
63
## for each component separately and is alligned with the components, so if the
 
64
## input's components were columns vectors the output will be a row vector and
 
65
## vice versa.
 
66
##
 
67
## @strong{Algorithms}
 
68
##
 
69
## The algorithms for this functions have been taken from the TISEAN package.
 
70
## @end deftypefn
 
71
 
 
72
## Author: Piotr Held <pjheld@gmail.com>.
 
73
## This function is based on timerev of TISEAN 3.0.1
 
74
## https://github.com/heggus/Tisean"
 
75
 
 
76
function output = timerev (S, delay);
 
77
 
 
78
  # Input validation
 
79
  if (nargin != 1 && nargin != 2)
 
80
    print_usage;
 
81
  endif
 
82
 
 
83
  if ((!ismatrix (S)) || (!isreal(S)))
 
84
    error ('Octave:invalid-input-arg', "S is not a realmatrix");
 
85
  endif
 
86
 
 
87
  # If no delay was given assign 1
 
88
  if (nargin == 1)
 
89
    delay = 1;
 
90
  endif
 
91
 
 
92
  # Verify delay is a positive integer
 
93
  isPositiveInteger = @(x) isreal(x) && isscalar (x) && (x > 0) ...
 
94
                           && (x-round(x) == 0);
 
95
  if (!isPositiveInteger(delay))
 
96
    error ("Octave:invalid-input-arg", "delay must be a positive integer");
 
97
  endif
 
98
 
 
99
  # Verify the input series 'S' is at least as long as 'delay + 1'
 
100
  if (delay >= length (S))
 
101
    error ("Octave:invalid-input-arg", ["the length of the input series must ",...
 
102
                                        "be at least as long as 'delay + 1'"]);
 
103
  endif
 
104
 
 
105
  # Correct S to always have more rows than columns
 
106
  trnspsd = false;
 
107
  if (rows (S) < columns (S))
 
108
    S = S.';
 
109
    trnspsd = true;
 
110
  endif
 
111
 
 
112
  # Calculate output
 
113
  idx    = delay+1:rows(S);
 
114
  t2     = sum ((S(idx,:) - S(idx-delay,:)).^2);
 
115
  t3     = sum ((S(idx,:) - S(idx-delay,:)).^3);
 
116
  output = t3./t2;
 
117
 
 
118
  # Transpose output if input components were row vectors and not column vectors
 
119
  if (trnspsd)
 
120
    output = output.';
 
121
  endif
 
122
 
 
123
endfunction
 
124
 
 
125
%% Test output against TISEAN program 'timerev'
 
126
%!assert (timerev (henon (1000)(:,1),4),0.313235879,-1e-6)
 
127
%!assert (timerev (henon (1000), 4), [0.313235879, 0.312556326], -1e-6)
 
128
%!assert (timerev (henon(1000).', 4), [0.313235879; 0.312556326], -1e-6)
 
129
 
 
130
%% Testing input validation
 
131
%!error <positive> timerev (1, 0);
 
132
%!error <long> timerev ([1 2], 2);
 
133
 
 
134
%% Test if default values load properly
 
135
%!assert (timerev (henon (100)), timerev (henon (100),1))