-
Notifications
You must be signed in to change notification settings - Fork 26
/
makeFittingStruct_GLMbi.m
52 lines (48 loc) · 1.76 KB
/
makeFittingStruct_GLMbi.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
46
47
48
49
50
51
52
function gg = makeFittingStruct_GLMbi(krank,varargin)
% gg = makeFittingStruct_GLM(krank,dtStim,dtSp,klength,nkbasis,k0,nhbasis,lasthpeak)
%
% Initialize parameter struct for fitting generalized linear model (GLM),
% with bilinear (i.e., low-rank) parametrization stimulus filter
%
% Inputs:
% krank = rank of stim filter
% dtStim = bin size of stimulus (in s)
% dtSp = bin size for spike train (in s)
% klength = temporal length of stimulus filter, in # of bins (optional)
% nkbasis = # temporal basis vectors for stim filter (optional)
% nhbasis = # temporal basis vectors for spike hist filter h (optional)
% lasthpeak = time of peak of last h basis vector, in s (optional)
% k0 = initial estimate of stimulus filter (optional)
%
% Outputs:
% gg = GLM-fitting param structure
% ====================================================================
% Set up structure
gg = makeFittingStruct_GLM(varargin{:});
% Set additional fields needed by bilinear GLM
gg.kx = [];
gg.kxbas=[];
gg.kbasprs = [];
gg.krank = krank;
gg.ktype = 'bilinear';
% if initial filter passed in, use svd to set up initial K params (bilinearly parametrized)
if (length(varargin)>4) && (~isempty(varargin{3}))
k0 = varargin{5};
[u,s,v] = svd(k0);
mux = mean(k0)';
nkt = size(gg.ktbas,2);
nkx = size(k0,2);
kt = zeros(nkt,krank);
kx = zeros(nkx,krank);
for j = 1:krank;
if v(:,j)'*mux < 0 % Flip sign if shape looks opposite to k0
u(:,j) = -u(:,j); v(:,j) = -v(:,j);
end
kt(:,j) = (gg.ktbas'*gg.ktbas)\(gg.ktbas'*(u(:,j)*s(j,j)));
kx(:,j) = v(:,j);
end
gg.kxbas = speye(nkx);
gg.kt = kt;
gg.kx = kx;
gg.k = (gg.ktbas*gg.kt)*(gg.kxbas*gg.kx)';
end