-
Notifications
You must be signed in to change notification settings - Fork 56
/
gsp_regression_tv.m
140 lines (103 loc) · 3.14 KB
/
gsp_regression_tv.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
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
function sol = gsp_regression_tv(G, M, y, tau, param )
%GSP_REGRESSION_TV Regression using graph and TV
% Usage: sol = gsp_regression_tv(G, M, y, tau );
% sol = gsp_regression_tv(G, M, y, tau, param );
%
% Input parameters:
% G : Graph
% M : Mask (to determine which label is known)
% y : label (total size of the problem)
% tau : regularization parameter (weight for tv)
% param : optional structure of parameters
%
% Output parameters:
% sol : Solution of the problem
%
% This function solve the following problem
%
% .. argmin_x || M x - y ||_2^2 + tau || x ||_{G TV}
%
% If tau is set to zero, then the following problem is solved
%
% .. argmin_x || x ||_{G TV} s. t. M x - y = 0
%
% This function uses the UNLocBoX.
%
% See also: gsp_regression_tik gsp_classification_tv
% Author: Nathanael Perraudin
% Date : 24 July 2015
%% Optional parameters
if nargin < 5
param = struct;
end
if nargin < 4
tau = 0;
end
if ~isfield(param,'verbose'), param.verbose = 1; end
%% prepare the graph
G = gsp_estimate_lmax(G);
G = gsp_adj2vec(G);
paramsolver = param;
if tau > 0
M_op =@(x) bsxfun(@times, M, x);
fg.grad = @(x) 2 * M_op(M_op(x) - y);
fg.eval = @(x) norm(M_op(x) - y)^2;
fg.beta = 2;
% setting the function ftv
% paramtv.verbose = param.verbose-1;
% ftv.prox = @(x,T) gsp_prox_tv(x,tau*T,G,paramtv);
% ftv.eval = @(x) tau *sum(gsp_norm_tv(G,x));
paramtv.verbose = param.verbose-1;
ftv.prox = @(x,T) prox_l1(x,tau*T,paramtv);
ftv.eval = @(x) tau *sum(gsp_norm_tv(G,x));
ftv.L = @(x ) G.Diff * x;
ftv.Lt = @(x ) G.Diff' * x;
ftv.norm_L = G.lmax;
%% solve the problem
% setting the timestep
sol = solvep(y,{ftv,fg},paramsolver);
else
% % param_b2.verbose = param.verbose -1;
% % param_b2.y = y;
% % param_b2.A = @(x) M.*x;
% % param_b2.At = @(x) M.*x;
% % param_b2.tight = param.tight;
% % param_b2.epsilon = 0;
% % fproj.prox = @(x,T) proj_b2(x,T,param_b2);
% % fproj.eval = @(x) eps;
%
% fproj.prox = @(x,T) x - M.*x + M.*y;
% fproj.eval = @(x) eps;
%
% % setting the function ftv
%
% % paramtv.verbose = param.verbose-1;
% % ftv.prox = @(x,T) gsp_prox_tv(x,T,G,paramtv);
% % ftv.eval = @(x) sum(gsp_norm_tv(G,x));
%
% paramtv.verbose = param.verbose-1;
% ftv.prox = @(x,T) prox_l1(x,T,paramtv);
% ftv.eval = @(x) sum(gsp_norm_tv(G,x));
% ftv.L = @(x ) G.Diff * x;
% ftv.Lt = @(x ) G.Diff' * x;
% ftv.norm_L = G.lmax;
%
% %% solve the problem
% sol = solvep(y,{ftv,fproj},paramsolver);
A = G.Diff(:,~logical(M));
b = - G.Diff(:,logical(M)) * y(logical(M),:);
paramtv.verbose = param.verbose-1;
paramtv.y = b;
ftv.prox = @(x,T) prox_l1(x,T,paramtv);
ftv.eval = @(x) sum(sum(abs(A*x-b)));
ftv.L = @(x ) A * x;
ftv.Lt = @(x ) A' * x;
ftv.norm_L = G.lmax;
%% solve the problem
solt = solvep(y(~logical(M),:),{ftv},paramsolver);
sol = y;
sol(~logical(M),:) = solt;
end
% sol = sol(logical(1-M));
% sol = reshape(sol,[],size(M,2));
end