SHOGUN  4.0.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
general_altra.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 <string.h>
21 
22 void general_altra(double *x, double *v, int n, double *G, double *ind, int nodes, double mult)
23 {
24 
25  int i, j;
26  double lambda,twoNorm, ratio;
27 
28  /*
29  * test whether the first node is special
30  */
31  if ((int) ind[0]==-1){
32 
33  /*
34  *Recheck whether ind[1] equals to zero
35  */
36  if ((int) ind[1]!=-1){
37  printf("\n Error! \n Check ind");
38  exit(1);
39  }
40 
41  lambda=mult*ind[2];
42 
43  for(j=0;j<n;j++){
44  if (v[j]>lambda)
45  x[j]=v[j]-lambda;
46  else
47  if (v[j]<-lambda)
48  x[j]=v[j]+lambda;
49  else
50  x[j]=0;
51  }
52 
53  i=1;
54  }
55  else{
56  memcpy(x, v, sizeof(double) * n);
57  i=0;
58  }
59 
60  /*
61  * sequentially process each node
62  *
63  */
64  for(;i < nodes; i++){
65  /*
66  * compute the L2 norm of this group
67  */
68  twoNorm=0;
69  for(j=(int) ind[3*i]-1;j< (int) ind[3*i+1];j++)
70  twoNorm += x[(int) G[j]-1 ] * x[(int) G[j]-1 ];
71  twoNorm=sqrt(twoNorm);
72 
73  lambda=mult*ind[3*i+2];
74  if (twoNorm>lambda){
75  ratio=(twoNorm-lambda)/twoNorm;
76 
77  /*
78  * shrinkage this group by ratio
79  */
80  for(j=(int) ind[3*i]-1;j<(int) ind[3*i+1];j++)
81  x[(int) G[j]-1 ]*=ratio;
82  }
83  else{
84  /*
85  * threshold this group to zero
86  */
87  for(j=(int) ind[3*i]-1;j<(int) ind[3*i+1];j++)
88  x[(int) G[j]-1 ]=0;
89  }
90  }
91 }
92 
93 void general_altra_mt(double *X, double *V, int n, int k, double *G, double *ind, int nodes, double mult)
94 {
95  int i, j;
96 
97  double *x=(double *)malloc(sizeof(double)*k);
98  double *v=(double *)malloc(sizeof(double)*k);
99 
100  for (i=0;i<n;i++){
101  /*
102  * copy a row of V to v
103  *
104  */
105  for(j=0;j<k;j++)
106  v[j]=V[j*n + i];
107 
108  general_altra(x, v, k, G, ind, nodes, mult);
109 
110  /*
111  * copy the solution to X
112  */
113  for(j=0;j<k;j++)
114  X[j*n+i]=x[j];
115  }
116 
117  free(x);
118  free(v);
119 }
120 
121 void general_computeLambda2Max(double *lambda2_max, double *x, int n, double *G, double *ind, int nodes)
122 {
123  int i, j;
124  double twoNorm;
125 
126  *lambda2_max=0;
127 
128 
129 
130  for(i=0;i < nodes; i++){
131  /*
132  * compute the L2 norm of this group
133  */
134  twoNorm=0;
135  for(j=(int) ind[3*i]-1;j< (int) ind[3*i+1];j++)
136  twoNorm += x[(int) G[j]-1 ] * x[(int) G[j]-1 ];
137  twoNorm=sqrt(twoNorm);
138 
139  twoNorm=twoNorm/ind[3*i+2];
140 
141  if (twoNorm >*lambda2_max )
142  *lambda2_max=twoNorm;
143  }
144 }
145 
146 double general_treeNorm(double *x, int ldx, int n, double *G, double *ind, int nodes)
147 {
148 
149  int i, j;
150  double twoNorm, lambda;
151 
152  double tree_norm=0;
153 
154  /*
155  * test whether the first node is special
156  */
157  if ((int) ind[0]==-1){
158 
159  /*
160  *Recheck whether ind[1] equals to zero
161  */
162  if ((int) ind[1]!=-1){
163  printf("\n Error! \n Check ind");
164  exit(1);
165  }
166 
167  lambda=ind[2];
168 
169  for(j=0;j<n;j+=ldx){
170  tree_norm+=fabs(x[j]);
171  }
172 
173  tree_norm=tree_norm * lambda;
174 
175  i=1;
176  }
177  else{
178  i=0;
179  }
180 
181  /*
182  * sequentially process each node
183  *
184  */
185  for(;i < nodes; i++){
186  /*
187  * compute the L2 norm of this group
188 
189 */
190  twoNorm=0;
191  for(j=(int) ind[3*i]-1;j< (int) ind[3*i+1];j++)
192  twoNorm += x[(int) G[j]-1 ] * x[(int) G[j]-1 ];
193  twoNorm=sqrt(twoNorm);
194 
195  lambda=ind[3*i+2];
196 
197  tree_norm=tree_norm + lambda*twoNorm;
198  }
199  return tree_norm;
200 }
201 
202 double general_findLambdaMax(double *v, int n, double *G, double *ind, int nodes)
203 {
204 
205  int i;
206  double lambda=0,squaredWeight=0, lambda1,lambda2;
207  double *x=(double *)malloc(sizeof(double)*n);
208  double *ind2=(double *)malloc(sizeof(double)*nodes*3);
209  int num=0;
210 
211  for(i=0;i<n;i++){
212  lambda+=v[i]*v[i];
213  }
214 
215  if ( (int)ind[0]==-1 )
216  squaredWeight=n*ind[2]*ind[2];
217  else
218  squaredWeight=ind[2]*ind[2];
219 
220  for (i=1;i<nodes;i++){
221  squaredWeight+=ind[3*i+2]*ind[3*i+2];
222  }
223 
224  /* set lambda to an initial guess
225  */
226  lambda=sqrt(lambda/squaredWeight);
227 
228  /*
229  printf("\n\n lambda=%2.5f",lambda);
230  */
231 
232  /*
233  *copy ind to ind2,
234  *and scale the weight 3*i+2
235  */
236  for(i=0;i<nodes;i++){
237  ind2[3*i]=ind[3*i];
238  ind2[3*i+1]=ind[3*i+1];
239  ind2[3*i+2]=ind[3*i+2]*lambda;
240  }
241 
242  /* test whether the solution is zero or not
243  */
244  general_altra(x, v, n, G, ind2, nodes);
245  for(i=0;i<n;i++){
246  if (x[i]!=0)
247  break;
248  }
249 
250  if (i>=n) {
251  /*x is a zero vector*/
252  lambda2=lambda;
253  lambda1=lambda;
254 
255  num=0;
256 
257  while(1){
258  num++;
259 
260  lambda2=lambda;
261  lambda1=lambda1/2;
262  /* update ind2
263  */
264  for(i=0;i<nodes;i++){
265  ind2[3*i+2]=ind[3*i+2]*lambda1;
266  }
267 
268  /* compute and test whether x is zero
269  */
270  general_altra(x, v, n, G, ind2, nodes);
271  for(i=0;i<n;i++){
272  if (x[i]!=0)
273  break;
274  }
275 
276  if (i<n){
277  break;
278  /*x is not zero
279  *we have found lambda1
280  */
281  }
282  }
283 
284  }
285  else{
286  /*x is a non-zero vector*/
287  lambda2=lambda;
288  lambda1=lambda;
289 
290  num=0;
291  while(1){
292  num++;
293 
294  lambda1=lambda2;
295  lambda2=lambda2*2;
296  /* update ind2
297  */
298  for(i=0;i<nodes;i++){
299  ind2[3*i+2]=ind[3*i+2]*lambda2;
300  }
301 
302  /* compute and test whether x is zero
303  */
304  general_altra(x, v, n, G, ind2, nodes);
305  for(i=0;i<n;i++){
306  if (x[i]!=0)
307  break;
308  }
309 
310  if (i>=n){
311  break;
312  /*x is a zero vector
313  *we have found lambda2
314  */
315  }
316  }
317  }
318 
319  /*
320  printf("\n num=%d, lambda1=%2.5f, lambda2=%2.5f",num, lambda1,lambda2);
321  */
322 
323  while ( fabs(lambda2-lambda1) > lambda2 * 1e-10 ){
324 
325  num++;
326 
327  lambda=(lambda1+lambda2)/2;
328 
329  /* update ind2
330  */
331  for(i=0;i<nodes;i++){
332  ind2[3*i+2]=ind[3*i+2]*lambda;
333  }
334 
335  /* compute and test whether x is zero
336  */
337  general_altra(x, v, n, G, ind2, nodes);
338  for(i=0;i<n;i++){
339  if (x[i]!=0)
340  break;
341  }
342 
343  if (i>=n){
344  lambda2=lambda;
345  }
346  else{
347  lambda1=lambda;
348  }
349 
350  /*
351  printf("\n lambda1=%2.5f, lambda2=%2.5f",lambda1,lambda2);
352  */
353  }
354 
355  /*
356  printf("\n num=%d",num);
357 
358  printf(" lambda1=%2.5f, lambda2=%2.5f",lambda1,lambda2);
359  */
360 
361  free(x);
362  free(ind2);
363 
364  return lambda2;
365 }
366 
367 double general_findLambdaMax_mt(double *V, int n, int k, double *G, double *ind, int nodes)
368 {
369  int i, j;
370 
371  double *v=(double *)malloc(sizeof(double)*k);
372  double lambda;
373 
374  double lambdaMax=0;
375 
376  for (i=0;i<n;i++){
377  /*
378  * copy a row of V to v
379  *
380  */
381  for(j=0;j<k;j++)
382  v[j]=V[j*n + i];
383 
384  lambda = general_findLambdaMax(v, k, G, ind, nodes);
385 
386  /*
387  printf("\n lambda=%5.2f",lambda);
388  */
389 
390 
391  if (lambda>lambdaMax)
392  lambdaMax=lambda;
393  }
394 
395  /*
396  printf("\n *lambdaMax=%5.2f",*lambdaMax);
397  */
398 
399  free(v);
400  return lambdaMax;
401 }

SHOGUN Machine Learning Toolbox - Documentation