SHOGUN  3.2.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
flsa.cpp
Go to the documentation of this file.
1 /* This program is free software: you can redistribute it and/or modify
2  * it under the terms of the GNU General Public License as published by
3  * the Free Software Foundation, either version 3 of the License, or
4  * (at your option) any later version.
5  *
6  * This program is distributed in the hope that it will be useful,
7  * but WITHOUT ANY WARRANTY; without even the implied warranty of
8  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
9  * GNU General Public License for more details.
10  *
11  * You should have received a copy of the GNU General Public License
12  * along with this program. If not, see <http://www.gnu.org/licenses/>.
13  *
14  * Copyright (C) 2009 - 2012 Jun Liu and Jieping Ye
15  */
16 
19 #include <stdlib.h>
20 #include <stdio.h>
21 #include <time.h>
22 #include <math.h>
23 
24 void flsa(double *x, double *z, double *infor,
25  double * v, double *z0,
26  double lambda1, double lambda2, int n,
27  int maxStep, double tol, int tau, int flag)
28 {
29 
30  int i, nn=n-1, m;
31  double zMax, temp;
32  double *Av, *g, *s;
33  int iterStep, numS;
34  double gap;
35  double *zz = NULL; /*to replace z0, so that z0 shall not revised after */
36 
37 
38  Av=(double *) malloc(sizeof(double)*nn);
39 
40  /*
41  Compute Av= A*v (n=4, nn=3)
42  A= [ -1 1 0 0;
43  0 -1 1 0;
44  0 0 -1 1]
45  */
46 
47  for (i=0;i<nn; i++)
48  Av[i]=v[i+1]-v[i];
49 
50  /*
51  Sovlve the linear system via Thomas's algorithm or Rose's algorithm
52  B * z0 = Av
53  */
54 
55  Thomas(&zMax, z, Av, nn);
56 
57  /*
58  We consider two cases:
59  1) lambda2 >= zMax, which leads to a solution with same entry values
60  2) lambda2 < zMax, which needs to first run sfa, and then perform soft thresholding
61  */
62 
63  /*
64  First case: lambda2 >= zMax
65  */
66  if (lambda2 >= zMax){
67 
68  temp=0;
69  m=n%5;
70  if (m!=0){
71  for (i=0;i<m;i++)
72  temp+=v[i];
73  }
74  for (i=m;i<n;i+=5){
75  temp += v[i] + v[i+1] + v[i+2] + v[i+3] + v[i+4];
76  }
77  temp/=n;
78  /* temp is the mean value of v*/
79 
80 
81  /*
82  soft thresholding by lambda1
83  */
84  if (temp> lambda1)
85  temp= temp-lambda1;
86  else
87  if (temp < -lambda1)
88  temp= temp+lambda1;
89  else
90  temp=0;
91 
92  m=n%7;
93  if (m!=0){
94  for (i=0;i<m;i++)
95  x[i]=temp;
96  }
97  for (i=m;i<n;i+=7){
98  x[i] =temp;
99  x[i+1] =temp;
100  x[i+2] =temp;
101  x[i+3] =temp;
102  x[i+4] =temp;
103  x[i+5] =temp;
104  x[i+6] =temp;
105  }
106 
107  gap=0;
108 
109  free(Av);
110 
111  if (infor)
112  {
113  infor[0]= gap;
114  infor[1]= 0;
115  infor[2]=zMax;
116  infor[3]=0;
117  }
118 
119  return;
120  }
121 
122 
123  /*
124  Second case: lambda2 < zMax
125 
126  We need to call sfa for computing x, and then do soft thresholding
127 
128  Before calling sfa, we need to allocate memory for g and s,
129  and initialize z and z0.
130  */
131 
132 
133  /*
134  Allocate memory for g and s
135  */
136 
137  g =(double *) malloc(sizeof(double)*nn),
138  s =(double *) malloc(sizeof(double)*nn);
139 
140 
141 
142  m=flag /10;
143  /*
144 
145  If m=0, then this shows that, z0 is a "good" starting point. (m=1-6)
146 
147  Otherwise (m=11-16), we shall set z as either the solution to the linear system.
148  or the zero point
149 
150 */
151  if (m==0){
152  for (i=0;i<nn;i++){
153  if (z0[i] > lambda2)
154  z[i]=lambda2;
155  else
156  if (z0[i]<-lambda2)
157  z[i]=-lambda2;
158  else
159  z[i]=z0[i];
160  }
161  }
162  else{
163  if (lambda2 >= 0.5 * zMax){
164  for (i=0;i<nn;i++){
165  if (z[i] > lambda2)
166  z[i]=lambda2;
167  else
168  if (z[i]<-lambda2)
169  z[i]=-lambda2;
170  }
171  }
172  else{
173  for (i=0;i<nn;i++)
174  z[i]=0;
175 
176  }
177  }
178 
179  flag=flag %10; /*
180  flag is now in [1:6]
181 
182  for sfa, i.e., flag in [1:4], we need initialize z0 with zero
183  */
184 
185  if (flag>=1 && flag<=4){
186  zz =(double *) malloc(sizeof(double)*nn);
187 
188  for (i=0;i<nn;i++)
189  zz[i]=0;
190  }
191 
192  /*
193  call sfa, sfa_one, or sfa_special to compute z, for finding the subgradient
194  and x
195  */
196 
197  if (flag==6)
198  iterStep=sfa_one(x, &gap, &numS,
199  z, v, Av,
200  lambda2, nn, maxStep,
201  s, g,
202  tol, tau);
203  else
204  if (flag==5)
205  iterStep=sfa_special(x, &gap, &numS,
206  z, v, Av,
207  lambda2, nn, maxStep,
208  s, g,
209  tol, tau);
210  else{
211  iterStep=sfa(x, &gap, &numS,
212  z, zz, v, Av,
213  lambda2, nn, maxStep,
214  s, g,
215  tol,tau, flag);
216 
217  free (zz);
218  /*free the variable zz*/
219  }
220 
221 
222  /*
223  soft thresholding by lambda1
224  */
225 
226  for(i=0;i<n;i++)
227  if (x[i] > lambda1)
228  x[i]-=lambda1;
229  else
230  if (x[i]<-lambda1)
231  x[i]+=lambda1;
232  else
233  x[i]=0;
234 
235 
236  free(Av);
237  free(g);
238  free(s);
239 
240  if (infor)
241  {
242  infor[0]=gap;
243  infor[1]=iterStep;
244  infor[2]=zMax;
245  infor[3]=numS;
246  }
247 }
248 

SHOGUN Machine Learning Toolbox - Documentation