v0.3.5:
[mldemos:mldemos.git] / _AlgorithmsPlugins / Maximizers / maximizeParticles.cpp
1 /*********************************************************************\r
2 MLDemos: A User-Friendly visualization toolkit for machine learning\r
3 Copyright (C) 2010  Basilio Noris\r
4 Contact: mldemos@b4silio.com\r
5 \r
6 This library is free software; you can redistribute it and/or\r
7 modify it under the terms of the GNU Lesser General Public\r
8 License as published by the Free Software Foundation; either\r
9 version 2.1 of the License, or (at your option) any later version.\r
10 \r
11 This library is distributed in the hope that it will be useful,\r
12 but WITHOUT ANY WARRANTY; without even the implied warranty of\r
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU\r
14 Library General Public License for more details.\r
15 \r
16 You should have received a copy of the GNU Lesser General Public\r
17 License along with this library; if not, write to the Free\r
18 Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.\r
19 *********************************************************************/\r
20 #include "public.h"\r
21 #include "basicMath.h"\r
22 #include "maximizeParticles.h"\r
23 #include <QDebug>\r
24 \r
25 using namespace std;\r
26 \r
27 MaximizeParticles::MaximizeParticles()\r
28         : particleCount(20)\r
29 {\r
30         dim = 2;\r
31         maximum.resize(dim);\r
32         FOR(d,dim) maximum[d] = rand()/(float)RAND_MAX;\r
33 }\r
34 \r
35 MaximizeParticles::~MaximizeParticles()\r
36 {\r
37         KILL(data);\r
38 }\r
39 \r
40 void MaximizeParticles::SetParams(int particleCount, float variance, bool bAdaptive)\r
41 {\r
42         this->particleCount = particleCount;\r
43         this->variance = variance;\r
44         this->bAdaptive = bAdaptive;\r
45         this->particleCount = particleCount;\r
46 }\r
47 \r
48 void MaximizeParticles::Draw(QPainter &painter)\r
49 {\r
50         painter.setPen(QPen(Qt::black, 1.5));\r
51         painter.setBrush(Qt::NoBrush);\r
52         FOR(i, visited.size())\r
53         {\r
54                 QPointF point(visited[i][0]*w, visited[i][1]*h);\r
55                 painter.drawEllipse(point, 3, 3);\r
56         }\r
57 \r
58         painter.setPen(QPen(Qt::black, 1.5));\r
59         FOR(i, history.size()-1 )\r
60         {\r
61                 QPointF point(history[i][0]*w, history[i][1]*h);\r
62                 QPointF pointNext(history[i+1][0]*w, history[i+1][1]*h);\r
63                 painter.setBrush(Qt::NoBrush);\r
64                 painter.drawLine(point, pointNext);\r
65                 painter.setBrush(Qt::white);\r
66                 painter.drawEllipse(point, 4, 4);\r
67         }\r
68 \r
69         // draw the current particles\r
70         FOR(i, particles.size())\r
71         {\r
72                 fvec sample = particles[i];\r
73                 QPointF point(sample[0]*w, sample[1]*h);\r
74                 int radius = 2 + weights[i]*5;\r
75                 painter.setBrush(Qt::green);\r
76                 painter.drawEllipse(point, radius, radius);\r
77         }\r
78         // we draw the current maximum\r
79         QPointF point(history[history.size()-1][0]*w, history[history.size()-1][1]*h);\r
80         painter.setBrush(QColor(255*(1-historyValue[history.size()-1]), 255, 255*(1-historyValue[history.size()-1])));\r
81         painter.drawEllipse(point, 5, 5);\r
82 }\r
83 \r
84 void MaximizeParticles::Train(float *dataMap, fVec size, fvec startingPoint)\r
85 {\r
86         w = size.x;\r
87         h = size.y;\r
88         if(data) delete [] data;\r
89         data = new float[w*h];\r
90         memcpy(data, dataMap, w*h*sizeof(float));\r
91         bConverged = false;\r
92         if(startingPoint.size())\r
93         {\r
94                 maximum = startingPoint;\r
95                 int xIndex = startingPoint[0]*w;\r
96                 int yIndex = startingPoint[1]*h;\r
97                 int index = min(w*h, max(0, yIndex*w + xIndex));\r
98                 float value = data[index];\r
99                 maximumValue = value;\r
100                 history.push_back(maximum);\r
101                 historyValue.push_back(value);\r
102                 //qDebug() << "Starting maximization at " << maximum[0] << " " << maximum[1];\r
103         }\r
104         particles.clear();\r
105         weights.clear();\r
106         fvec sample; sample.resize(dim);\r
107         FOR(i, particleCount)\r
108         {\r
109                 FOR(d, dim) sample[d] = drand48();\r
110                 particles.push_back(sample);\r
111                 weights.push_back(1.f/particleCount);\r
112         }\r
113 }\r
114 \r
115 fvec MaximizeParticles::Test( const fvec &sample)\r
116 {\r
117         if(bConverged) return maximum;\r
118 \r
119         float decay = 0.2f;\r
120         float totalWeights= 0;\r
121         FOR(i, particleCount)\r
122         {\r
123                 // first we guess the next pose for each particle\r
124                 particles[i] += RandN(dim, 0, variance*variance);\r
125                 // we compute the weights\r
126                 weights[i] = weights[i] *(1-decay) + GetValue(particles[i])*decay;\r
127                 totalWeights += weights[i];\r
128                 //visited.push_back(particles[i]);\r
129         }\r
130 \r
131         // we sort the particles by weight\r
132         vector< pair<double, u32> > fits;\r
133         FOR(i, weights.size()) fits.push_back(pair<double, u32>(weights[i], i));\r
134         sort(fits.begin(), fits.end(), greater< pair<double,u32> >());\r
135         vector<fvec> newParticles;\r
136         fvec newWeights;\r
137         FOR(i, fits.size())\r
138         {\r
139                 newParticles.push_back(particles[fits[i].second]);\r
140                 newWeights.push_back(weights[fits[i].second]);\r
141         }\r
142         particles = newParticles;\r
143         weights = newWeights;\r
144 \r
145         // we compute the result\r
146         maximum = particles[0];\r
147         maximumValue = weights[0];\r
148 \r
149         if(bAdaptive)\r
150         {\r
151                 float decay = 0.1;\r
152                 variance = variance*(1-decay) + (1 - maximumValue*maximumValue)*0.1*decay;\r
153         }\r
154 \r
155 /*\r
156         maximum.clear();\r
157         maximum.resize(dim,0);\r
158         FOR(i, particleCount)\r
159         {\r
160                 FOR(d, dim) maximum[d] = weights[i] * particles[i][d];\r
161         }\r
162 */\r
163 \r
164         history.push_back(maximum);\r
165         historyValue.push_back(maximumValue);\r
166 \r
167         // we normalize probabilities by weight value\r
168         double weightSum = 0;\r
169         FOR(i, weights.size())  weightSum += weights[i];\r
170         dvec probs;\r
171         double weightCounter = 0;\r
172         FOR(i, weights.size())\r
173         {\r
174                 weightCounter += weights[i]/weightSum;\r
175                 probs.push_back(weightCounter);\r
176         }\r
177 \r
178         // we resample the particles\r
179         newParticles.clear();\r
180         newWeights.clear();\r
181         FOR(i, particleCount)\r
182         {\r
183                 double r = drand48();\r
184                 u32 j=0;\r
185                 for (j=0; j<probs.size() && r > probs[j]; j++);\r
186                 u32 index = j<probs.size() ? j : 0;\r
187                 newParticles.push_back(particles[index]);\r
188                 newWeights.push_back(weights[index]);\r
189         }\r
190         particles = newParticles;\r
191         weights = newWeights;\r
192 \r
193         // if particles have gone wild we resample them\r
194         float degeneracyThreshold = 0.000000001;\r
195         FOR(i, particleCount)\r
196         {\r
197                 if( weights[i] < degeneracyThreshold )\r
198                 {\r
199                         FOR(d, dim) particles[i][d] = drand48();\r
200                         weights[i] = 1.f / particleCount;\r
201                 }\r
202         }\r
203         return maximum;\r
204 }\r
205 \r
206 fvec MaximizeParticles::Test(const fVec &sample)\r
207 {\r
208         return Test((fvec)sample);\r
209 }\r
210 \r
211 char *MaximizeParticles::GetInfoString()\r
212 {\r
213         char *text = new char[1024];\r
214         sprintf(text, "Genetic Algorithm\n");\r
215         return text;\r
216 }\r