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

[scipy-svn@scip...]

Author: cdavid
Date: 2007-07-02 04:49:12 -0500 (Mon, 02 Jul 2007)
New Revision: 3137

Modified:
   trunk/Lib/sandbox/pyem/densities.py
   trunk/Lib/sandbox/pyem/examples/pdfestimation1d.py
Log:
Clean up for densities.py

Modified: trunk/Lib/sandbox/pyem/densities.py
===================================================================
--- trunk/Lib/sandbox/pyem/densities.py	2007-07-02 09:31:18 UTC (rev 3136)
+++ trunk/Lib/sandbox/pyem/densities.py	2007-07-02 09:49:12 UTC (rev 3137)
@@ -1,7 +1,7 @@
 #! /usr/bin/python
 #
 # Copyrighted David Cournapeau
-# Last Change: Sat Jun 30 04:00 PM 2007 J
+# Last Change: Mon Jul 02 06:00 PM 2007 J
 """This module implements various basic functions related to multivariate
 gaussian, such as pdf estimation, confidence interval/ellipsoids, etc..."""
 
@@ -294,55 +294,19 @@
     mu = N.atleast_2d(mu)
     va = N.atleast_2d(va)
 
-    K = N.shape(mu)[0]
+    k = N.shape(mu)[0]
     n = N.shape(data)[0]
     d = N.shape(mu)[1]
     
-    y = N.zeros((K, n))
+    y = N.zeros((k, n))
     if N.size(mu) == N.size(va):
-        for i in range(K):
+        for i in range(k):
             y[i] = gauss_den(data, mu[i, :], va[i, :], log)
         return y.T
     else:
-        for i in range(K):
+        for i in range(k):
             y[i] = gauss_den(data, mu[i, :], va[d*i:d*i+d, :], log)
         return y.T
 
 if __name__ == "__main__":
     pass
-    ## import pylab
-
-    ## #=========================================
-    ## # Test plotting a simple diag 2d variance:
-    ## #=========================================
-    ## va  = N.array([5, 3])
-    ## mu  = N.array([2, 3])
-
-    ## # Generate a multivariate gaussian of mean mu and covariance va
-    ## X       = randn(1e3, 2)
-    ## Yc      = N.dot(N.diag(N.sqrt(va)), X.transpose())
-    ## Yc      = Yc.transpose() + mu
-
-    ## # Plotting
-    ## Xe, Ye  = gauss_ell(mu, va, npoints = 100)
-    ## pylab.figure()
-    ## pylab.plot(Yc[:, 0], Yc[:, 1], '.')
-    ## pylab.plot(Xe, Ye, 'r')
-
-    ## #=========================================
-    ## # Test plotting a simple full 2d variance:
-    ## #=========================================
-    ## va  = N.array([[0.2, 0.1],[0.1, 0.5]])
-    ## mu  = N.array([0, 3])
-
-    ## # Generate a multivariate gaussian of mean mu and covariance va
-    ## X       = randn(1e3, 2)
-    ## Yc      = N.dot(lin.cholesky(va), X.transpose())
-    ## Yc      = Yc.transpose() + mu
-
-    ## # Plotting
-    ## Xe, Ye  = gauss_ell(mu, va, npoints = 100, level=0.95)
-    ## pylab.figure()
-    ## pylab.plot(Yc[:, 0], Yc[:, 1], '.')
-    ## pylab.plot(Xe, Ye, 'r')
-    ## pylab.show()

Modified: trunk/Lib/sandbox/pyem/examples/pdfestimation1d.py
===================================================================
--- trunk/Lib/sandbox/pyem/examples/pdfestimation1d.py	2007-07-02 09:31:18 UTC (rev 3136)
+++ trunk/Lib/sandbox/pyem/examples/pdfestimation1d.py	2007-07-02 09:49:12 UTC (rev 3137)
@@ -51,11 +51,17 @@
 
 mbic = N.argmax(bc)
 
-# Below is code to display a figure with histogram and best model (in the BIC sense)
-# pdf, with the BIC as a function of the number of components on the right.
+# Below is code to display a figure with histogram and best model (in the BIC
+# sense) pdf, with the BIC as a function of the number of components on the
+# right.
 P.figure(figsize = [12, 7])
+#---------------------------
+# histogram + pdf estimation
+#---------------------------
 P.subplot(1, 2, 1)
 h = gml[mbic].plot1d(gpdf=True)
+# You can manipulate the differents parts of the plot through the returned
+# handles
 h['gpdf'][0].set_linestyle('-')
 h['gpdf'][0].set_label('pdf of the mixture')
 h['pdf'][0].set_label('pdf of individual component')
@@ -68,6 +74,9 @@
 P.hist(dt, 25, normed = 1, fill = False)
 P.xlabel('waiting time between consecutive eruption (in min)')
 
+#------------------------------------------
+# BIC as a function of number of components
+#------------------------------------------
 P.subplot(1, 2, 2)
 P.plot(N.arange(1, 8), bc, 'o:')
 P.xlabel("number of components")