[Scipy-svn] r2275 - in trunk/Lib/sandbox/pyem: . pyem

scipy-svn at scipy.org scipy-svn at scipy.org
Thu Oct 12 08:47:06 CDT 2006


Author: cdavid
Date: 2006-10-12 08:46:56 -0500 (Thu, 12 Oct 2006)
New Revision: 2275

Modified:
   trunk/Lib/sandbox/pyem/Changelog
   trunk/Lib/sandbox/pyem/pyem/__init__.py
   trunk/Lib/sandbox/pyem/pyem/densities.py
   trunk/Lib/sandbox/pyem/pyem/gauss_mix.py
   trunk/Lib/sandbox/pyem/pyem/gmm_em.py
   trunk/Lib/sandbox/pyem/pyem/kmean.py
   trunk/Lib/sandbox/pyem/pyem/profile_densities.py
   trunk/Lib/sandbox/pyem/pyem/profile_gmm.py
Log:
[pyem @ david at ar.media.kyoto-u.ac.jp-20060817094851-03f3aebbbf6b4161]
Adapt to numpy 1.0b3SVN
David Cournapeau <david at ar.media.kyoto-u.ac.jp> | 2006-08-17 18:48:51 +0900 (Thu, 17 Aug 2006)

Modified: trunk/Lib/sandbox/pyem/Changelog
===================================================================
--- trunk/Lib/sandbox/pyem/Changelog	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/Changelog	2006-10-12 13:46:56 UTC (rev 2275)
@@ -1,3 +1,11 @@
+pyem (0.5.1) Thu, 17 Aug 2006 11:54:41 +0900
+
+	* put version to 0.5.1
+	* Update to last numpy (change axis args between 1.0b1 and 1.0b2,
+	change type args of ones and zeros)
+
+-- David Cournapeau <david at ar.media.kyoto-u.ac.jp> 
+
 pyem (0.5) Fri, 04 Aug 2006 23:10:37 +0900
 
 	* put version to 0.5.0

Modified: trunk/Lib/sandbox/pyem/pyem/__init__.py
===================================================================
--- trunk/Lib/sandbox/pyem/pyem/__init__.py	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/pyem/__init__.py	2006-10-12 13:46:56 UTC (rev 2275)
@@ -1,7 +1,7 @@
 #! /usr/bin/env python
-# Last Change: Fri Aug 04 11:00 PM 2006 J
+# Last Change: Thu Aug 17 03:00 PM 2006 J
 
-version = '0.5.0'
+version = '0.5.1'
 
 from gauss_mix import GmParamError, GM
 from gmm_em import GmmParamError, GMM

Modified: trunk/Lib/sandbox/pyem/pyem/densities.py
===================================================================
--- trunk/Lib/sandbox/pyem/pyem/densities.py	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/pyem/densities.py	2006-10-12 13:46:56 UTC (rev 2275)
@@ -1,7 +1,7 @@
 #! /usr/bin/python
 #
 # Copyrighted David Cournapeau
-# Last Change: Fri Aug 04 11:00 PM 2006 J
+# Last Change: Thu Aug 17 03:00 PM 2006 J
 
 import numpy as N
 import numpy.linalg as lin
@@ -146,7 +146,7 @@
 
     # # Slow version
     # n       = N.size(x, 0)
-    # y       = N.zeros(n, float)
+    # y       = N.zeros(n)
     # for i in range(n):
     #     y[i] = N.dot(x[i,:],
     #              N.dot(inva, N.transpose(x[i,:])))
@@ -201,7 +201,7 @@
     mu  = mu[dim]
     if mode == 'diag':
         va      = va[dim]
-        elps    = N.outer(mu, N.ones(npoints, float))
+        elps    = N.outer(mu, N.ones(npoints))
         elps    += N.dot(N.diag(N.sqrt(va)), circle)
     elif mode == 'full':
         va  = va[c,:][:,c]
@@ -213,7 +213,7 @@
         #   - va = cova' * cova (matlab)
         # So take care when comparing results with matlab !
         cova    = lin.cholesky(va)
-        elps    = N.outer(mu, N.ones(npoints, float))
+        elps    = N.outer(mu, N.ones(npoints))
         elps    += N.dot(cova, circle)
     else:
         raise DenParam("var mode not recognized")

Modified: trunk/Lib/sandbox/pyem/pyem/gauss_mix.py
===================================================================
--- trunk/Lib/sandbox/pyem/pyem/gauss_mix.py	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/pyem/gauss_mix.py	2006-10-12 13:46:56 UTC (rev 2275)
@@ -1,5 +1,5 @@
 # /usr/bin/python
-# Last Change: Fri Aug 04 07:00 PM 2006 J
+# Last Change: Thu Aug 17 03:00 PM 2006 J
 
 # Module to implement GaussianMixture class.
 
@@ -54,12 +54,12 @@
 
         # Init to 0 all parameters, with the right dimensions.
         # Not sure this is useful in python from an efficiency POV ?
-        self.w   = N.zeros(k, float)
-        self.mu  = N.zeros((k, d), float)
+        self.w   = N.zeros(k)
+        self.mu  = N.zeros((k, d))
         if mode == 'diag':
-            self.va  = N.zeros((k, d), float)
+            self.va  = N.zeros((k, d))
         elif mode == 'full':
-            self.va  = N.zeros((k * d, d), float)
+            self.va  = N.zeros((k * d, d))
 
         self.is_valid   = False
 
@@ -96,7 +96,7 @@
             X   = self.mu[S, :]  + X * N.sqrt(self.va[S,:])
         elif self.mode == 'full':
             # Faster:
-            cho = N.zeros((self.k, self.va.shape[1], self.va.shape[1]), float)
+            cho = N.zeros((self.k, self.va.shape[1], self.va.shape[1]))
             for i in range(self.k):
                 # Using cholesky looks more stable than sqrtm; sqrtm is not
                 # available in numpy anyway, only in scipy...
@@ -110,14 +110,13 @@
 
         return X
 
-    def conf_ellipses(self, c = [0, 1], npoints = 100):
+    def conf_ellipses(self, *args, **kargs):
         """Returns a list of confidence ellipsoids describing the Gmm
-        defined by mu and va. c is the dimension we are projecting
-        the variances on a 2d space. For now, the confidence level
-        is fixed to 0.39.
-        
+        defined by mu and va. Check densities.gauss_ell for details
+
         Returns:
-            -Xe:    a list of x coordinates for the ellipses
+            -Xe:    a list of x coordinates for the ellipses (Xe[i] is
+            the array containing x coordinates of the ith Gaussian)
             -Ye:    a list of y coordinates for the ellipses
 
         Example:
@@ -134,22 +133,19 @@
             Will plot samples X draw from the mixture model, and
             plot the ellipses of equi-probability from the mean with
             fixed level of confidence 0.39.  """
-        # TODO: adjustable level (to do in gauss_ell). 
-        # For now, a level of 0.39 means that we draw
-        # ellipses for 1 standard deviation. 
         Xe  = []
         Ye  = []   
         if self.mode == 'diag':
             for i in range(self.k):
-                xe, ye  = densities.gauss_ell(self.mu[i,:], self.va[i,:], dim = c, 
-                        npoints = npoints)
+                xe, ye  = densities.gauss_ell(self.mu[i,:], self.va[i,:], 
+                        *args, **kargs)
                 Xe.append(xe)
                 Ye.append(ye)
         elif self.mode == 'full':
             for i in range(self.k):
                 xe, ye  = densities.gauss_ell(self.mu[i,:], 
-                        self.va[i*self.d:i*self.d+self.d,:], dim = c, 
-                        npoints = npoints)
+                        self.va[i*self.d:i*self.d+self.d,:], 
+                        *args, **kargs)
                 Xe.append(xe)
                 Ye.append(ye)
 
@@ -165,7 +161,7 @@
         Returns: w, mu, va
         """
         w   = abs(randn(nc))
-        w   = w / sum(w)
+        w   = w / sum(w, 0)
 
         mu  = spread * randn(nc, d)
         if varmode == 'diag':
@@ -222,7 +218,7 @@
     """
         
     # Check that w is valid
-    if N.fabs(N.sum(w)  - 1) > MAX_DEV:
+    if N.fabs(N.sum(w, 0)  - 1) > MAX_DEV:
         raise GmmParamError('weight does not sum to 1')
     
     if not len(w.shape) == 1:
@@ -266,7 +262,7 @@
     #   - d = dimension
     #   - mode : mode of covariance matrices
     d       = 5
-    k       = 5
+    k       = 4
     mode    = 'full'
     nframes = 1e3
 
@@ -283,9 +279,10 @@
 
     P.plot(X[:, 0], X[:, 1], '.', label = '_nolegend_')
 
-    # Real confidence ellipses with level 0.39
-    Xre, Yre  = gm.conf_ellipses()
-    P.plot(Xre[0], Yre[0], 'g', label = 'true confidence ellipsoides')
+    # Real confidence ellipses with confidence level 
+    level       = 0.39
+    Xre, Yre    = gm.conf_ellipses(level = level)
+    P.plot(Xre[0], Yre[0], 'g', label = 'true confidence ellipsoides at level %.2f' % level)
     for i in range(1,k):
         P.plot(Xre[i], Yre[i], 'g', label = '_nolegend_')
 

Modified: trunk/Lib/sandbox/pyem/pyem/gmm_em.py
===================================================================
--- trunk/Lib/sandbox/pyem/pyem/gmm_em.py	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/pyem/gmm_em.py	2006-10-12 13:46:56 UTC (rev 2275)
@@ -1,5 +1,5 @@
 # /usr/bin/python
-# Last Change: Fri Aug 04 07:00 PM 2006 J
+# Last Change: Thu Aug 17 03:00 PM 2006 J
 
 import numpy as N
 import numpy.linalg as lin
@@ -60,15 +60,15 @@
 
         (code, label)   = kmean(data, init, niter)
 
-        w   = N.ones(k, float) / k
+        w   = N.ones(k) / k
         mu  = code.copy()
         if self.gm.mode == 'diag':
-            va = N.zeros((k, d), float)
+            va = N.zeros((k, d))
             for i in range(k):
                 for j in range(d):
                     va[i,j] = N.cov(data[N.where(label==i), j], rowvar = 0)
         elif self.gm.mode == 'full':
-            va  = N.zeros((k*d, d), float)
+            va  = N.zeros((k*d, d))
             for i in range(k):
                 va[i*d:i*d+d,:] = \
                     N.cov(data[N.where(label==i)], rowvar = 0)
@@ -84,7 +84,7 @@
         k   = self.gm.k
         d   = self.gm.d
         if mode == 'diag':
-            w   = N.ones(k, float) / k
+            w   = N.ones(k) / k
             mu  = randn(k, d)
             va  = N.fabs(randn(k, d))
         else:
@@ -142,11 +142,11 @@
         d       = self.gm.d
         n       = data.shape[0]
         invn    = 1.0/n
-        mGamma  = N.sum(gamma)
+        mGamma  = N.sum(gamma, axis = 0)
 
         if self.gm.mode == 'diag':
-            mu  = N.zeros((k, d), float)
-            va  = N.zeros((k, d), float)
+            mu  = N.zeros((k, d))
+            va  = N.zeros((k, d))
             gamma   = gamma.transpose()
             for c in range(k):
                 # x       = N.sum(N.outerproduct(gamma[:, k], 
@@ -161,18 +161,18 @@
             w   = invn * mGamma
 
         elif self.gm.mode == 'full':
-            mu  = N.zeros((k, d), float)
-            va  = N.zeros((k*d, d), float)
+            mu  = N.zeros((k, d))
+            va  = N.zeros((k*d, d))
 
             for c in range(k):
                 x   = N.sum(N.outerproduct(gamma[:, c], 
-                            N.ones((1, d), float)) * data)
-                xx  = N.zeros((d, d), float)
+                            N.ones((1, d))) * data, axis = 0)
+                xx  = N.zeros((d, d))
                 
                 # This should be much faster than recursing on n...
                 for i in range(d):
                     for j in range(d):
-                        xx[i,j] = N.sum(data[:,i] * data[:,j] * gamma[:,c])
+                        xx[i,j] = N.sum(data[:,i] * data[:,j] * gamma[:,c], axis = 0)
 
                 mu[c,:] = x / mGamma[c]
                 va[c*d:c*d+d,:] = xx  / mGamma[c] - \
@@ -211,12 +211,12 @@
         model.init(data)
 
         # Likelihood is kept
-        like    = N.zeros(niter, float)
+        like    = N.zeros(niter)
 
         # Em computation, with computation of the likelihood
         for i in range(niter):
             g, tgd      = model.sufficient_statistics(data)
-            like[i]     = N.sum(N.log(N.sum(tgd, 1)))
+            like[i]     = N.sum(N.log(N.sum(tgd, 1)), axis = 0)
             model.update_em(data, g)
 
         return like
@@ -232,7 +232,7 @@
     n   = data.shape[0]
     d   = data.shape[1]
     
-    y   = N.zeros((n, K), float)
+    y   = N.zeros((n, K))
     if mu.size == va.size:
         for i in range(K):
             y[:, i] = densities.gauss_den(data, mu[i, :], va[i, :])
@@ -288,12 +288,12 @@
 
     # The actual EM, with likelihood computation
     niter   = 10
-    like    = N.zeros(niter, float)
+    like    = N.zeros(niter)
 
     print "computing..."
     for i in range(niter):
         g, tgd  = gmm.sufficient_statistics(data)
-        like[i] = N.sum(N.log(N.sum(tgd, 1)))
+        like[i] = N.sum(N.log(N.sum(tgd, 1)), axis = 0)
         gmm.update_em(data, g)
     # # Alternative form, by using EM class: as the EM class
     # # is quite rudimentary now, it is not very useful, just save

Modified: trunk/Lib/sandbox/pyem/pyem/kmean.py
===================================================================
--- trunk/Lib/sandbox/pyem/pyem/kmean.py	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/pyem/kmean.py	2006-10-12 13:46:56 UTC (rev 2275)
@@ -1,5 +1,5 @@
 # /usr/bin/python
-# Last Change: Thu Jul 13 07:00 PM 2006 J
+# Last Change: Thu Aug 17 03:00 PM 2006 J
 
 import numpy as N
 
@@ -20,14 +20,14 @@
         msg = "data and init centers do not have same dimensions..."
         raise GmmParamError(msg)
     
-    code    = N.asarray(init.copy(), float)
+    code    = N.asarray(init.copy())
     for i in range(iter):
         # Compute the nearest neighbour for each obs
         # using the current code book
         label   = _vq(data, code)
         # Update the code by computing centroids using the new code book
         for j in range(k):
-            code[j,:] = N.mean(data[N.where(label==j)]) 
+            code[j,:] = N.mean(data[N.where(label==j)], axis=0) 
 
     return code, label
 

Modified: trunk/Lib/sandbox/pyem/pyem/profile_densities.py
===================================================================
--- trunk/Lib/sandbox/pyem/pyem/profile_densities.py	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/pyem/profile_densities.py	2006-10-12 13:46:56 UTC (rev 2275)
@@ -45,4 +45,3 @@
     import pstats
     p = pstats.Stats('gdenprof')
     print p.sort_stats('cumulative').print_stats(20)
-

Modified: trunk/Lib/sandbox/pyem/pyem/profile_gmm.py
===================================================================
--- trunk/Lib/sandbox/pyem/pyem/profile_gmm.py	2006-10-12 13:46:47 UTC (rev 2274)
+++ trunk/Lib/sandbox/pyem/pyem/profile_gmm.py	2006-10-12 13:46:56 UTC (rev 2275)
@@ -39,7 +39,7 @@
 
     # The actual EM, with likelihood computation
     niter   = 10
-    like    = N.zeros(niter, float)
+    like    = N.zeros(niter)
 
     print "computing..."
     for i in range(niter):



More information about the Scipy-svn mailing list