Home | Develop | Download | Contact
test_working_FittingGradientIGSoft.cpp
1 
20 #include <iostream>
21 #include <Pds/Ra>
22 #include <Pds/Ml>
23 #include <Pds/Sp>
24 
25 namespace Pds{
26 namespace LogisticModel{
27 
28 Pds::Vector FittingGradientIGSoft( Pds::IterationConf &Conf,
29  const Pds::Matrix &X,
30  const Pds::Vector &Y,
31  const Pds::Vector &W0);
32 
33 Pds::Vector GradientCostInformationGainSoft(const Pds::Vector &W,
34  const Pds::Matrix &X,
35  const Pds::Vector &Y,
36  double h);
37 
38 Pds::Vector GradientCostInformationGainSoft2( const Pds::Vector &W,
39  const Pds::Matrix &X,
40  const Pds::Vector &Y);
41 
42 double CostInformationGainSoft( const Pds::Vector &W,
43  const Pds::Matrix &X,
44  const Pds::Vector &Y);
45 
46 
47 }//LogisticModel
48 }//Pds
49 
50 
51 double Pds::LogisticModel::CostInformationGainSoft( const Pds::Vector &W,
52  const Pds::Matrix &X,
53  const Pds::Vector &Y)
54 {
55  if(W.IsEmpty()) return Pds::Ra::Nan;
56  if(X.IsEmpty()) return Pds::Ra::Nan;
57  if(Y.IsEmpty()) return Pds::Ra::Nan;
58 
59  if(X.Nlin()!=Y.Nlin()) return Pds::Ra::Nan;
60  if(X.Ncol()!=(W.Nlin()-1)) return Pds::Ra::Nan;
61 
62  Pds::Vector Yh=Y.Geq(0.5);
63 
64  Pds::Vector z=Pds::LogisticModel::Classify(W,X);
65 
66  unsigned int L=Yh.Nlin();
67 
68  double Nt=Yh.Sum();
69  double Na=z.Sum();
70  double N1=Yh.Dot(z);
71 
72  /*
73  // Information gain
74  double Et=Pds::Hb(Nt/L);
75  double E1=Pds::Hb(N1/Na);
76  double val=(Nt-N1)/(L-Na);
77  if(val>1.0) val=1.0;
78  if(val<0) val=0;
79  if(L==Na) val=1;
80  double E0=Pds::Hb(val);
81  */
82 
83  // Gini index
84  double q=2.0;
85  double F=Pds::qHbn(0.5,2.0);
86  double Et=Pds::qHbn(Nt/L,q)/F;
87  double E1=Pds::qHbn(N1/Na,q)/F;
88  double val=(Nt-N1)/(L-Na);
89  if(val>1.0) val=1.0;
90  if(val<0) val=0;
91  if(L==Na) val=1;
92  double E0=Pds::qHbn(val,q)/F;
93 
94 
95 
96  return Et-(Na/L)*E1-((L-Na)/L)*E0;
97 }
98 
99 
100 
101 Pds::Vector Pds::LogisticModel::GradientCostInformationGainSoft(const Pds::Vector &W,
102  const Pds::Matrix &X,
103  const Pds::Vector &Y,
104  double h)
105 {
106  if(W.IsEmpty()) return Pds::Vector();
107  if(X.IsEmpty()) return Pds::Vector();
108  if(Y.IsEmpty()) return Pds::Vector();
109 
110  if(X.Nlin()!=Y.Nlin()) return Pds::Vector();
111  if(X.Ncol()!=(W.Nlin()-1)) return Pds::Vector();
112 
113  if(h==0.0) return Pds::Vector();
114 
115  unsigned int N=W.Nel();
116 
117  Pds::Vector D(N);
118  Pds::Vector Whp(N);
119  Pds::Vector Whm(N);
120  double dn;
121 
122  Whp.Copy(W);
123  Whm.Copy(W);
124 
125  for(unsigned int n=0;n<N;n++)
126  {
127  Whp.SetRaw(n,W.GetRaw(n)+h);
128  Whm.SetRaw(n,W.GetRaw(n)-h);
129 
130  dn=( Pds::LogisticModel::CostInformationGainSoft(Whp,X,Y)
131  -Pds::LogisticModel::CostInformationGainSoft(Whm,X,Y) )/(2*h);
132 
133  D.SetRaw(n,dn);
134 
135  Whp.SetRaw(n,W.GetRaw(n));
136  Whm.SetRaw(n,W.GetRaw(n));
137  }
138  return D;
139 }
140 
141 Pds::Vector Pds::LogisticModel::GradientCostInformationGainSoft2(const Pds::Vector &W,
142  const Pds::Matrix &X,
143  const Pds::Vector &Y)
144 {
145  if(W.IsEmpty()) return Pds::Vector();
146  if(X.IsEmpty()) return Pds::Vector();
147  if(Y.IsEmpty()) return Pds::Vector();
148 
149  if(X.Nlin()!=Y.Nlin()) return Pds::Vector();
150  if(X.Ncol()!=(W.Nlin()-1)) return Pds::Vector();
151 
152  Pds::Vector Yt=Y.Geq(0.5);
153  unsigned int N=Yt.Sum();
154  unsigned int L=Yt.Nel();
155  Pds::Matrix R=Pds::RegressorMatrix(X);
156  Pds::Vector z=Pds::Sigmoid(R*W);
157  Pds::Vector Dz=z.Product(1.0-z);
158 
159  double Na=z.Sum();
160  double N1=Yt.Dot(z);
161  double p1=N1/Na;
162  double p2=(N-N1)/(L-Na);
163 
164  double h1=Pds::Hb(p1);
165  double h2=Pds::Hb(p2);
166  double Dh1=-Pds::Logit2(p1);
167  double Dh2=-Pds::Logit2(p2);
168  /*
169  double q=4.0;
170  double K=Pds::qHbn(0.5,q);
171  double h1=Pds::qHbn(p1,q)/K;
172  double h2=Pds::qHbn(p2,q)/K;
173  double Dh1=Pds::qDHbn(p1,q)/K;
174  double Dh2=Pds::qDHbn(p2,q)/K;
175  */
176  double factor1=(-h1+h2+p1*Dh1-p2*Dh2)/L;
177  double factor2=(-Dh1+Dh2)/L;
178  Pds::Vector DIG= factor1*R.TMul(Dz)+factor2*R.TMul(Dz.Product(Yt));
179 
180  std::cout<<"\n";
181  std::cout<<"N:"<<N<<"\n";
182  std::cout<<"L:"<<L<<"\n";
183  std::cout<<"N1:"<<N1<<"\t N1<Na True positive in Na side\n";
184  std::cout<<"Na:"<<Na<<"\t elems. classify as 1\n";
185  std::cout<<"p1:"<<p1<<"\n";
186  std::cout<<"p2:"<<p2<<"\n";
187  std::cout<<"factor1:"<<factor1<<"\n";
188  std::cout<<"factor2:"<<factor2<<"\n";
189  std::cout<<"IGs:"<<Pds::LogisticModel::CostInformationGainSoft(W,X,Y)<<"\n";
190  std::cout<<"IGh:"<<Pds::LogisticModel::CostInformationGain(W,X,Y)<<"\n";
191  DIG.T().Print("DIG: ");
192 
193  return DIG;
194 }
195 
196 
197 
198 void print_IterationConf_IG_soft_init_data(Pds::IterationConf Conf,std::string str="")
199 {
200  double Alpha = Conf.GetAlpha();
201  double MinError= Conf.GetMinError();
202  double MaxIter = Conf.GetMaxIter();
203 
204  std::ios_base::fmtflags f( std::cout.flags() );
205 
206  std::cout.precision(4);
207  std::cout<<"┌──────────────────────────────┐"<<std::endl;
208  std::cout<<"│ "<<std::left<<std::setw(22)<<str<<" │"<<std::endl;
209  std::cout<<"│ IterationConf init data │"<<std::endl;
210  std::cout<<"├──────────────────────────────┤"<<std::endl;
211  std::cout<<"│ Alpha: "<<std::setw(10)<<Alpha <<" │"<<std::endl;
212  std::cout<<"│ MinErrorΔ: "<<std::setw(10)<<std::scientific<<MinError<<std::defaultfloat<<" │"<<std::endl;
213  std::cout<<"│ MaxIter: "<<std::setw(10)<<MaxIter<<" │"<<std::endl;
214  std::cout<<"└──────────────────────────────┘"<<std::endl;
215 
216  std::cout.flags( f );
217 }
218 
219 void print_IterationConf_IG_soft_end_data(Pds::IterationConf Conf,std::string str="")
220 {
221  std::ios_base::fmtflags f( std::cout.flags() );
222 
223  double Alpha = Conf.GetAlpha();
224  std::cout.precision(3);
225  std::cout<<"┌──────────────────────────────┐"<<std::endl;
226  std::cout<<"│ "<<std::left<<std::setw(22)<<str<<" │"<<std::endl;
227  std::cout<<"│ IterationConf end data │"<<std::endl;
228  std::cout<<"├──────────────────────────────┤"<<std::endl;
229  std::cout<<"│ Alpha: "<<std::setw(10)<<Alpha<<" │"<<std::endl;
230  std::cout<<"│ ErrorΔ: "<<std::setw(10)<<std::scientific<<Conf.LastError<<std::defaultfloat<<" │"<<std::endl;
231  std::cout<<"│ LastIter: "<<std::setw(10)<<Conf.LastIter<<" │"<<std::endl;
232  std::cout<<"└──────────────────────────────┘"<<std::endl;
233 
234  std::cout.flags( f );
235 }
236 void print_iterative_IG_soft_data(unsigned int iter,double Alpha,double last_error,double cost_error,bool END)
237 {
238 
239  std::ios_base::fmtflags f( std::cout.flags() );
240 
241  std::cout.precision(3);
242  std::cout<<std::left
243  <<"┌──────────────────────────────┐"<<std::endl
244  <<"│ Iter: "<<std::setw(10)<<iter<<" │"<<std::endl
245  <<"│ Alpha: "<<std::setw(10)<<Alpha<<" │"<<std::endl
246  <<"│ ErrorΔ: "<<std::setw(10)<<std::scientific<<last_error<<" │"<<std::endl
247  <<"│ IG: "<<std::setw(10)<<cost_error<<std::defaultfloat<<" │"<<std::endl
248  <<"└──────────────────────────────┘"<<std::endl;
249  if(END==false) std::cout<<"\x1b[A\r"<<"\x1b[A\r"<<"\x1b[A\r"<<"\x1b[A\r"<<"\x1b[A\r"<<"\x1b[A\r";
250 
251  std::cout.flags( f );
252 }
253 
254 Pds::Vector Pds::LogisticModel::FittingGradientIGSoft( Pds::IterationConf &Conf,
255  const Pds::Matrix &X,
256  const Pds::Vector &Y,
257  const Pds::Vector &W0)
258 {
259  if(W0.IsEmpty()) return Pds::Matrix();
260  if(X.IsEmpty()) return Pds::Matrix();
261  if(Y.IsEmpty()) return Pds::Matrix();
262 
263  if(X.Nlin()!=Y.Nlin()) return Pds::Matrix();
264  if(X.Ncol()!=(W0.Nlin()-1)) return Pds::Matrix();
265 
266  unsigned int YSUM=Y.Geq(0.5).Sum();
267  if(YSUM==Y.Nel()) return Pds::LogisticModel::GetW0CornerMeanMethod(X);
268  if(YSUM==0) return -Pds::LogisticModel::GetW0CornerMeanMethod(X);
269 
270  double Alpha=Conf.GetAlpha();
271  //double Gamma=Conf.GetGamma();
272  //double Lambda=Conf.GetLambda();
273  double MinError=Conf.GetMinError();
274  double MaxIter=Conf.GetMaxIter();
275 
276  pds_print_error_message("Esta función aun no trabaja bien, es demasiado abrupta");
277 
278  if(Conf.Show) print_IterationConf_IG_soft_init_data(Conf,"Gradient IG");
279 
280  Pds::Vector W=W0;
281  Pds::Vector Wopt=W;
282  Pds::Vector DIG,DCE;
283  double Costmax=-1;
284  double Cost;
285 
286 
287  Pds::Fir FIR5(Pds::Vector(5,1.0)/5);
288  //double fir5;
289  //Pds::Fir FIR10(Pds::Vector(10,1.0)/10);
290  //double fir10;
291 
292  double last_error=1000, before_error=0,delta_error=1;;
293  unsigned int iter=0;
294 
295  Pds::Matrix R=Pds::RegressorMatrix(X);
296  Pds::Vector DY,DW;
297  Pds::Matrix I=Pds::Eye(W.Nel());I.SetRaw(0,0,0);
298  Pds::Vector Yt;
299 
300 
301  std::vector<double> ig;
302  std::vector<double> err;
303 
304  double A=W.Norm()*0.1+0.01;
305 
306  do{
307  //double h;
308  //h=W.RMS()*0.5;
309  //DIG=Pds::LogisticModel::GradientCostInformationGainSoft(W,X,Y,h);
310  DW=-Pds::LogisticModel::GradientCostInformationGainSoft2(W,X,Y);
311  DW.T().Print("DW: ");
312  DW=(A/DW.Norm())*DW;
313  DW.T().Print("DW: ");
314 
315  W=W-Alpha*DW;
316 
317  W.T().Print(" W: ");
318 
319  Pds::Octave::Plot::PointsX2DY(X,Pds::Sigmoid(R*W),"testando.m","test_"+std::to_string(iter)+".png");
320 
321  Cost=Pds::LogisticModel::CostInformationGain(W,X,Y);
322 
323  before_error=last_error;
324 
325  last_error = 1.0-Cost;
326  //last_error += (Lambda*0.5/W.Nel())*(W.Dot(W)-W[0]*W[0]);
327 
328  delta_error=fabs(FIR5.Evaluate(last_error)-before_error);
329 
330  if(Cost>Costmax)
331  {
332  Costmax=Cost;
333  Wopt=W;
334  }
335 
336  ig.push_back(Cost);
337  err.push_back(last_error);
338 
339  iter++;
340 
341  //if(Conf.Show) print_iterative_IG_soft_data(iter,Alpha,delta_error,Cost,true);
342 
343  }while( (delta_error>MinError)&&(iter<MaxIter) );
344 
345  if(Conf.Show) print_iterative_IG_soft_data(iter,Alpha,delta_error,Cost,true);
346 
347  Conf.SetAlpha(Alpha);
348  Conf.LastError=delta_error;
349  Conf.LastIter=iter;
350 
351  Pds::Ra::SaveStdVector(ig,"ig.txt");
352  Pds::Ra::SaveStdVector(err,"err.txt");
353 
354  if(Conf.Show) print_IterationConf_IG_soft_end_data(Conf,"Gradient IG");
355 
356  double alpha=1.0;
357  /*
358  double s2;
359  Yt=R*Wopt;
360  double Delta=0.001;
361  Pds::Matrix Yo=Y.Geq(0.5)*(1.0-2*Delta)+Delta;
362  s2=Yt.SumSquare();
363  if(s2>0) alpha=Yt.Dot(Pds::Logit(Yo))/s2;
364  else alpha=1;
365  */
366  return Wopt*alpha;
367 }
368 
369 double average_distance_between_samples(const Pds::Matrix &X)
370 {
371  if(X.IsEmpty()) return Pds::Ra::Nan;
372 
373  Pds::Vector Min=X.MinInCols();
374  Pds::Vector Max=X.MaxInCols();
375 
376  unsigned int N=X.Ncol();
377  unsigned int L=X.Nlin();
378  unsigned int n;
379 
380  double S=1;
381  double d;
382  for(n=0;n<N;n++)
383  {
384  d=fabs(Max[n]-Min[n]);
385  if(d!=0) S=S*d;
386  }
387 
388  return pow(S/L,1.0/N);
389 }
390 
391 #include "extra_newhb.h"
392 
393 // ./test/test_working_FittingGradientIGSoft
394 int main(void)
395 {
396  Pds::IterationConf Conf;
397  Conf.Show=true;
398  Conf.SetMaxIter(100);
399  Conf.SetMinError(1e-07);
400  Conf.SetAlpha(0.5);
401  Conf.SetLambda(0.0);
402  Conf.SetGamma(0.1);
403  Pds::Ra::Randomize();
404 
405  // Generating data
406  /* // Data 1
407  Pds::Matrix X(Pds::Ra::TextFormat,"../test/hard-data-x.txt");
408  Pds::Vector Y(Pds::Ra::TextFormat,"../test/hard-data-y.txt");
409  X=X-X.MeanInCols();
410  */
411  // Data 2
412  unsigned int L=1000; Pds::Matrix X; Pds::Vector Y;
413  Pds::DataSet::Classification::LoadDataBand(L,X,Y);
414 
415  Pds::Octave::XLabel="x_1";
416  Pds::Octave::YLabel="x_2";
417  Pds::Octave::Plot::PointsX2DY(X,Y,"testando.m","test_working_FittingGradientIGSoft.png");
418 
419  /*
420  Pds::Vector W,DW;
421  for(unsigned int n=0;n<10;n++)
422  {
423  W=Pds::LogisticModel::GetW0MeanMethod(X);
424  double h=W.RMS()*0.001;
425  DW=Pds::LogisticModel::GradientCostInformationGainSoft(W,X,Y,h);
426  DW.Normalize();
427  DW.T().Print("\nDW1: ");
428  DW=Pds::LogisticModel::GradientCostInformationGainSoft2(W,X,Y);
429  DW.Normalize();
430  DW.T().Print("DW2: ");
431  //W=Pds::LogisticModel::FittingGradientIGSoft(Conf,X,Y,W);
432  }
433  */
434 
435  /*
436  Pds::Vector W;
437  W=10*Pds::LogisticModel::GetW0MeanMethod(X);
438  W=Pds::LogisticModel::FittingGradientIGSoft(Conf,X,Y,W);
439  */
440  unsigned int N=128;
441  double b=4;
442  Pds::Vector w1=Pds::LinSpace(-b,b,N);
443  Pds::Vector w2=Pds::LinSpace(-b,b,N);
444  Pds::Matrix W1,W2;
445  Pds::Matrix IG(N,N);
446  Pds::Matrix IGsoft(N,N);
447  double val;
448 
449  std::cout<<"avr: "<<average_distance_between_samples(X)<<"\n";
450 
451  Pds::Meshgrid (w1,w2,W1,W2);
452 
453  for(unsigned int i=0;i<N;i++)
454  for(unsigned int j=0;j<N;j++)
455  {
456  Pds::Vector w({0.5,w1[i],w2[j]});
457  val=Pds::LogisticModel::CostInformationGain( w,X,Y);
458  IG.SetRaw(i,j,val);
459  val=Pds::LogisticModel::CostXqEntropy( w,X,Y,q_factor_val); //w*(100/w.Norm())
460  IGsoft.SetRaw(i,j,val);
461  }
462 
463  // No treinar bias
464  // Agrega ||W-w_i|| e inicia desde best ortogonal
465 
466  Pds::Octave::Plot::SurfC(W1,W2,(1-IG),"testando.m","ig_matrix.png");
467  Pds::Octave::Plot::SurfC(W1,W2,IGsoft,"testando_soft.m","ig_soft_matrix.png");
468 
469  return 0;
470 }
Pds::IterationConf
La clase tipo Pds::IterationConf . Esta clase genera una matriz de Nlin lineas y 1 columna....
Definition: IterationConf.hpp:57
Pds::IterationConf::SetLambda
bool SetLambda(double Lambda)
Coloca el valor lambda.
Pds::IterationConf::SetAlpha
bool SetAlpha(double Alpha)
Coloca el valor alpha.
Pds::IterationConf::GetMinError
double GetMinError(void) const
Devuelve el valor MinError.
Pds::IterationConf::GetAlpha
double GetAlpha(void) const
Devuelve el valor alpha.
Pds::IterationConf::SetMaxIter
bool SetMaxIter(unsigned int MaxIter)
Coloca el valor MaxIter.
Pds::IterationConf::SetMinError
bool SetMinError(double MinError)
Coloca el valor MinError.
Pds::IterationConf::SetGamma
bool SetGamma(double Lambda)
Coloca el valor gamma.
Pds::IterationConf::GetMaxIter
unsigned int GetMaxIter(void) const
Devuelve el valor MaxIter.
Pds::DataSet::Classification::LoadDataBand
void LoadDataBand(unsigned int L, Pds::Matrix &X, Pds::Vector &Y)
Clasificacion de datos separados por mas de una curva.
Pds::LogisticModel::GetW0CornerMeanMethod
Pds::Vector GetW0CornerMeanMethod(const Pds::Matrix &X, double Delta=0.001)
Obtiene de forma rapida un vector inicial para usar en regresion logistica.
Pds::LogisticModel::CostInformationGain
double CostInformationGain(const Pds::Vector &W, const Pds::Matrix &X, const Pds::Vector &Y)
Calculo de costo.
Pds::LogisticModel::Classify
Pds::Vector Classify(const Pds::Vector &W, const Pds::Matrix &X)
Calculo del resultado del clasificador.
Pds::LogisticModel::CostInformationGainSoft
double CostInformationGainSoft(const Pds::Vector &W, const Pds::Matrix &X, const Pds::Vector &Y)
Definition: test_working_FittingGradientIGSoft.cpp:51
Pds::LogisticModel::FittingGradientIGSoft
Pds::Vector FittingGradientIGSoft(Pds::IterationConf &Conf, const Pds::Matrix &X, const Pds::Vector &Y, const Pds::Vector &W0)
Definition: test_working_FittingGradientIGSoft.cpp:254
Pds::LogisticModel::GradientCostInformationGainSoft2
Pds::Vector GradientCostInformationGainSoft2(const Pds::Vector &W, const Pds::Matrix &X, const Pds::Vector &Y)
Definition: test_working_FittingGradientIGSoft.cpp:141
Pds::LogisticModel::CostXqEntropy
double CostXqEntropy(const Pds::Vector &W, const Pds::Matrix &X, const Pds::Vector &Y, double q)
Definition: extra_newhb.h:32
Pds::LogisticModel::GradientCostInformationGainSoft
Pds::Vector GradientCostInformationGainSoft(const Pds::Vector &W, const Pds::Matrix &X, const Pds::Vector &Y, double h)
Definition: test_working_FittingGradientIGSoft.cpp:101
Pds
Nombre de espacio para Pds (Procesamiento Digital de Senales)
Definition: ClassificationMetrics.hpp:42

Enlaces de interés

HomePage Bazaar Download Bug report Ayuda Developer Feed
Manual de referencia generado el Domingo, 28 de Noviembre de 2021 13:03:56 para la biblioteca libpdsmlmm usando  Generado el Domingo, 28 de Noviembre de 2021 13:03:56 para libpdsmlmm por   doxygen 1.9.2