-
Notifications
You must be signed in to change notification settings - Fork 10
/
Copy pathiss_ds.m
45 lines (38 loc) · 1.31 KB
/
iss_ds.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
% function from State-Space Granger Causality Matlab® Toolbox
% based on "Barnett, L. and Seth, A. K.,Granger causality for state-space models, 91(4),2015"
% downloaded from http://users.sussex.ac.uk/~lionelb/
function [Ak,Kk,Vk] = iss_ds(A,C,K,V,k)
% Calculate innovations form state space parameters for a downsampled state
% space model
%
% A,C,K,V - innovations form state space parameters
% k - downsample factor (positive integer)
%
% Ak,Kk,Vk - innovations form state space parameters for downsampled model
%
% See V. Solo, arXiv:1501.04663v1, Sec. 4.
[m,m1] = size(A); assert(m1 == m);
[n,m1] = size(C); assert(m1 == m);
[m1,n1] = size(K); assert(n1 == n && m1 == m);
[n1,n2] = size(V); assert(n1 == n && n2 == n);
assert(isscalar(k) && isnumeric(k) && k == floor(k) && k > 0);
if k == 1 % nothing to do
Ak = A;
Kk = K;
Vk = V;
return
end
Ak = eye(m); % Ak = A^0
AkKSQRTV = K*chol(V,'lower');
Qk = AkKSQRTV*AkKSQRTV';
for p = 1:k-1
Ak = Ak*A; % Ak = A^p
AkKSQRTV = Ak*AkKSQRTV;
Qk = Qk + AkKSQRTV*AkKSQRTV';
end
Sk = Ak*K*V; % Ak = A^(k-1)
Ak = Ak*A; % Ak = A^k
% We now have general form state space parameters for the downsampled model
% (note: Ck = C and Rk = V for any k). We solve the DARE to get remaining
% innovations form parameters.
[Kk,Vk] = iss_ss2iss(Ak,C,Qk,V,Sk);