Completed
Push — 0.7.dev ( bb7559...c599bb )
by Andrei
57s
created

ema.__log_likelihood()   A

Complexity

Conditions 3

Size

Total Lines 11

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 3
dl 0
loc 11
rs 9.4285
c 0
b 0
f 0
1
"""!
2
3
@brief Cluster analysis algorithm: Expectation-Maximization Algorithm (EMA).
4
@details Implementation based on article:
5
         - 
6
7
@authors Andrei Novikov ([email protected])
8
@date 2014-2017
9
@copyright GNU Public License
10
11
@cond GNU_PUBLIC_LICENSE
12
    PyClustering is free software: you can redistribute it and/or modify
13
    it under the terms of the GNU General Public License as published by
14
    the Free Software Foundation, either version 3 of the License, or
15
    (at your option) any later version.
16
    
17
    PyClustering is distributed in the hope that it will be useful,
18
    but WITHOUT ANY WARRANTY; without even the implied warranty of
19
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
20
    GNU General Public License for more details.
21
    
22
    You should have received a copy of the GNU General Public License
23
    along with this program.  If not, see <http://www.gnu.org/licenses/>.
24
@endcond
25
26
"""
27
28
29
import numpy;
0 ignored issues
show
Configuration introduced by
The import numpy could not be resolved.

This can be caused by one of the following:

1. Missing Dependencies

This error could indicate a configuration issue of Pylint. Make sure that your libraries are available by adding the necessary commands.

# .scrutinizer.yml
before_commands:
    - sudo pip install abc # Python2
    - sudo pip3 install abc # Python3
Tip: We are currently not using virtualenv to run pylint, when installing your modules make sure to use the command for the correct version.

2. Missing __init__.py files

This error could also result from missing __init__.py files in your module folders. Make sure that you place one file in each sub-folder.

Loading history...
30
31
from pyclustering.cluster import cluster_visualizer;
32
from pyclustering.cluster.center_initializer import kmeans_plusplus_initializer;
33
from pyclustering.cluster.kmeans import kmeans;
34
35
from pyclustering.utils import pi, calculate_ellipse_description;
36
37
import matplotlib.pyplot as plt;
0 ignored issues
show
Configuration introduced by
The import matplotlib.pyplot could not be resolved.

This can be caused by one of the following:

1. Missing Dependencies

This error could indicate a configuration issue of Pylint. Make sure that your libraries are available by adding the necessary commands.

# .scrutinizer.yml
before_commands:
    - sudo pip install abc # Python2
    - sudo pip3 install abc # Python3
Tip: We are currently not using virtualenv to run pylint, when installing your modules make sure to use the command for the correct version.

2. Missing __init__.py files

This error could also result from missing __init__.py files in your module folders. Make sure that you place one file in each sub-folder.

Loading history...
38
from matplotlib import patches;
0 ignored issues
show
Configuration introduced by
The import matplotlib could not be resolved.

This can be caused by one of the following:

1. Missing Dependencies

This error could indicate a configuration issue of Pylint. Make sure that your libraries are available by adding the necessary commands.

# .scrutinizer.yml
before_commands:
    - sudo pip install abc # Python2
    - sudo pip3 install abc # Python3
Tip: We are currently not using virtualenv to run pylint, when installing your modules make sure to use the command for the correct version.

2. Missing __init__.py files

This error could also result from missing __init__.py files in your module folders. Make sure that you place one file in each sub-folder.

Loading history...
39
40
41
42
def gaussian(data, mean, covariance):
43
    dimension = float(len(data[0]));
44
 
45
    if (dimension != 1.0):
46
        inv_variance = numpy.linalg.pinv(covariance);
47
    else:
48
        inv_variance = 1.0 / covariance;
49
    
50
    divider = (pi * 2.0) ** (dimension / 2.0) * numpy.sqrt(numpy.linalg.norm(covariance));
51
    right_const = 1.0 / divider;
52
     
53
    result = [];
54
     
55
    for point in data:
56
        mean_delta = point - mean;
57
        point_gaussian = right_const * numpy.exp( -0.5 * mean_delta.dot(inv_variance).dot(numpy.transpose(mean_delta)) );
58
        result.append(point_gaussian);
59
     
60
    return result;
61
62
63
64
class ema_observer:
65
    def __init__(self):
66
        self.__means_evolution = [];
67
        self.__covariances_evolution = [];
68
        self.__clusters_evolution = [];
69
70
71
    def get_iterations(self):
72
        return len(self.__means);
0 ignored issues
show
Bug introduced by
The Instance of ema_observer does not seem to have a member named __means.

This check looks for calls to members that are non-existent. These calls will fail.

The member could have been renamed or removed.

Loading history...
73
74
75
    def get_means(self):
76
        return self.__means_evolution;
77
78
79
    def get_covariances(self):
80
        return self.__covariances_evolution;
81
82
83
    def notify(self, means, covariances, clusters):
84
        self.__means_evolution.append(means);
85
        self.__covariances_evolution.append(covariances);
86
        self.__clusters_evolution.append(clusters);
87
88
89
90
class ema_visualizer:
91
    @staticmethod
92
    def show_clusters(clusters, sample, covariances, means, display = True):
93
        visualizer = cluster_visualizer();
94
        visualizer.append_clusters(clusters, sample);
95
        figure = visualizer.show(display = False);
96
        
97
        if (len(sample[0]) == 2):
98
            ema_visualizer.__draw_ellipses(figure, visualizer, clusters, covariances, means);
99
100
        if (display is True): 
101
            plt.show();
102
103
        return figure;
104
105
106
    @staticmethod
107
    def __draw_ellipses(figure, visualizer, clusters, covariances, means):
108
        print(len(clusters));
109
        print([len(cluster) for cluster in clusters]);
110
        print(clusters);
111
        
112
        ax = figure.get_axes()[0];
113
        
114
        for index in range(len(clusters)):
115
            angle, width, height = calculate_ellipse_description(covariances[index]);
116
            color = visualizer.get_cluster_color(index, 0);
117
            
118
            ema_visualizer.__draw_ellipse(ax, means[index][0], means[index][1], angle, width, height, color);
119
120
121
    @staticmethod
122
    def __draw_ellipse(ax, x, y, angle, width, height, color):
123
        ellipse = patches.Ellipse((x, y), width, height, alpha=0.2, angle=angle, linewidth=2, fill=True, zorder=2, color=color);
124
        ax.add_patch(ellipse);
125
126
127
class ema:
128
    def __init__(self, data, amount_clusters, means = None, variances = None, observer = None, tolerance = 0.00001):
129
        self.__data = numpy.array(data);
130
        self.__amount_clusters = amount_clusters;
131
        self.__tolerance = tolerance;
132
        self.__observer = observer;
133
        
134
        self.__means = means;
135
        if ((means is None) or (variances is None)):
136
            self.__means, self.__variances = self.__get_initial_parameters(data, amount_clusters);
137
        
138
        self.__rc = [ [0.0] * len(self.__data) for _ in range(amount_clusters) ];
139
        self.__pic = [1.0] * amount_clusters;
140
        self.__clusters = [];
141
        self.__gaussians = [ [] for _ in range(amount_clusters) ];
142
        self.__stop = False;
143
144
145
    def process(self):
146
        self.__clusters = None;
147
        
148
        previous_likelihood = -200000;
149
        current_likelihood = -100000;
150
        
151
        while( (self.__stop is False) and (abs(previous_likelihood - current_likelihood) > self.__tolerance) ):
152
            self.__expectation_step();
153
            self.__maximization_step();
154
            
155
            previous_likelihood = current_likelihood;
156
            current_likelihood = self.__log_likelihood();
157
            self.__stop = self.__get_stop_condition();
158
        
159
        self.__clusters = self.__extract_clusters();
160
161
162
    def get_clusters(self):
163
        return self.__clusters;
164
165
166
    def get_centers(self):
167
        return self.__means;
168
169
170
    def get_covariances(self):
171
        return self.__variances;
172
173
174
    def __notify(self):
175
        if (self.__observer is not None):
176
            clusters = self.__extract_clusters();
177
            self.__notify(self.__means, self.__variances, clusters);
0 ignored issues
show
Bug introduced by
There seem to be too many positional arguments for this method call.
Loading history...
178
179
180
    def __extract_clusters(self):
181
        clusters = [ [] for _ in range(self.__amount_clusters) ];
182
        for index_point in range(len(self.__data)):
183
            candidates = [];
184
            for index_cluster in range(self.__amount_clusters):
185
                candidates.append((index_cluster, self.__rc[index_cluster][index_point]));
186
            
187
            index_winner = max(candidates, key = lambda candidate : candidate[1])[0];
188
            clusters[index_winner].append(index_point);
189
        
190
        clusters = [ cluster for cluster in clusters if len(cluster) > 0 ];
191
        return clusters;
192
193
194
    def __log_likelihood(self):
195
        likelihood = 0.0;
196
        
197
        for index_point in range(len(self.__data)):
198
            particle = 0.0;
199
            for index_cluster in range(self.__amount_clusters):
200
                particle += self.__pic[index_cluster] * self.__gaussians[index_cluster][index_point];
201
            
202
            likelihood += numpy.log(particle);
203
        
204
        return likelihood;
205
206
207
    def __probabilities(self, index_cluster, index_point):
208
        divider = 0.0;
209
        for i in range(self.__amount_clusters):
210
            divider += self.__pic[i] * self.__gaussians[i][index_point];
211
        
212
        rc = self.__pic[index_cluster] * self.__gaussians[index_cluster][index_point] / divider;
213
        return rc;
214
215
216
    def __expectation_step(self):
217
        for index in range(self.__amount_clusters):
218
            self.__gaussians[index] = gaussian(self.__data, self.__means[index], self.__variances[index]);
219
        
220
        for index_cluster in range(self.__amount_clusters):
221
            for index_point in range(len(self.__data)):
222
                self.__rc[index_cluster][index_point] = self.__probabilities(index_cluster, index_point);
223
224
225
    def __maximization_step(self):
226
        self.__pic = [];
227
        self.__means = [];
228
        self.__variances = [];
229
        
230
        amount_impossible_clusters = 0;
231
        
232
        for index_cluster in range(self.__amount_clusters):
233
            mc = numpy.sum(self.__rc[index_cluster]);
234
            
235
            if (mc == 0.0):
236
                amount_impossible_clusters += 1;
237
                continue;
238
            
239
            self.__pic.append( mc / len(self.__data) );
240
            self.__means.append( self.__update_mean(self.__rc[index_cluster], mc) );
241
            self.__variances.append( self.__update_covariance(self.__means[-1], self.__rc[index_cluster], mc) );
242
        
243
        self.__amount_clusters -= amount_impossible_clusters;
244
245
246
    def __get_stop_condition(self):
247
        for covariance in self.__variances:
248
            if (numpy.linalg.norm(covariance) == 0.0):
249
                return True;
250
        
251
        return False;
252
253
254
    def __update_covariance(self, means, rc, mc):
255
        covariance = 0.0;
256
        for index_point in range(len(self.__data)):
257
            deviation = numpy.array( [ self.__data[index_point] - means ]);
258
            covariance += rc[index_point] * deviation.T.dot(deviation);
259
        
260
        covariance = covariance / mc;
261
        return covariance;
262
263
264
    def __update_mean(self, rc, mc):
265
        mean = 0.0;
266
        for index_point in range(len(self.__data)):
267
            mean += rc[index_point] * self.__data[index_point];
268
        
269
        mean = mean / mc;
270
        return mean;
271
272
273
    def __get_initial_parameters(self, sample, amount_clusters):
274
        initial_centers = kmeans_plusplus_initializer(sample, amount_clusters).initialize();
275
        kmeans_instance = kmeans(sample, initial_centers, ccore = True);
276
        kmeans_instance.process();
277
        
278
        means = kmeans_instance.get_centers();
279
        
280
        covariances = [];
281
        initial_clusters = kmeans_instance.get_clusters();
282
        for initial_cluster in initial_clusters:
283
            cluster_sample = [ sample[index_point] for index_point in initial_cluster ];
284
            covariances.append(numpy.cov(cluster_sample, rowvar = False));
285
        
286
        return means, covariances;