[Neurodebian-devel] MDP patches in AFNI

Yaroslav Halchenko debian at onerussian.com
Tue Nov 20 17:41:44 UTC 2012


Thanks guys for the discussion.  Although I have not touched ICAs for a
while, having two separate limits for the coarse/fine optimizations
makes sense in general.  And it could be achieved without affecting
existing users with only 1 small change -- not relying only on
comparison of g and g_fine but rather on explicit setting for the
complementary (coarser) limit.  Since existing users had no clue about
Initial and Secondary, and messages "Secondary convergence..." without
ever seeing "Primary..." one, I suggest to make "primary_limit" a
"coarse_limit" and set it to None, by default disabling this
functionaliry.  Then whoever needs to enable this coarse_limit (e.g. for
meica) could do that.

I have pushed original patch + my changes on top resulting in listed
patch to:
https://github.com/yarikoptic/mdp-toolkit
and sent a PR

I then have altered fastica tests to sweep also through having some
coarse_limit set -- nothing crashes, although now it nearly doubles test run
time, so you might like to improve on that ? ;)  and I haven't tested with
python3 yet.

If this patch, in this or another shape gets accepted upstream
in MDP then we could update MDP's debian package and use system-wide MDP for
meica to avoid maintaining multiple copies/clones.

Cheers,

$> git diff origin/master..              
diff --git a/mdp/nodes/ica_nodes.py b/mdp/nodes/ica_nodes.py
index 6ba1a0a..31f1657 100644
--- a/mdp/nodes/ica_nodes.py
+++ b/mdp/nodes/ica_nodes.py
@@ -320,13 +320,14 @@ class FastICANode(ICANode):
     - 25.5.2005 now independent from scipy. Requires Numeric or numarray
     - 26.6.2006 converted to numpy
     - 14.9.2007 updated to Matlab version 2.5
+    - 26.6.2012 added ability to run two stages of optimization [PK]
     """
 
     def __init__(self, approach = 'defl', g = 'pow3', guess = None,
                  fine_g = 'pow3', mu = 1,
                  sample_size = 1, fine_tanh = 1, fine_gaus = 1,
-                 max_it = 1000, max_it_fine = 100,
-                 failures = 5, limit = 0.001, verbose = False,
+                 max_it = 5000, max_it_fine = 100,
+                 failures = 5, coarse_limit=None, limit = 0.001,  verbose = False,
                  whitened = False, white_comp = None, white_parm = None,
                  input_dim = None, dtype=None):
         """
@@ -346,6 +347,12 @@ def __init__(self, approach = 'defl', g = 'pow3', guess = None,
                       It is passed directly to the WhiteningNode constructor.
                       Ex: white_parm = { 'svd' : True }
 
+        coarse_limit -- initial convergence threshold, to switch to
+                        fine_g function (i.e. linear to non-linear)
+                        even before reaching the limit and final
+                        tuning. Set it to a value higher than limit to
+                        be in effect. PK 26-6-12.
+
         limit -- convergence threshold.
 
         Specific for FastICA:
@@ -417,6 +424,7 @@ def __init__(self, approach = 'defl', g = 'pow3', guess = None,
         self.fine_gaus = fine_gaus
         self.max_it = max_it
         self.max_it_fine = max_it_fine
+        self.coarse_limit = coarse_limit
         self.failures = failures
         self.guess = guess
 
@@ -458,6 +466,7 @@ def core(self, data):
                 guess = mult(guess, self.white.get_recmatrix(transposed=1))
 
         limit = self.limit
+        coarse_limit = self.coarse_limit
         max_it = self.max_it
         max_it_fine = self.max_it_fine
         failures = self.failures
@@ -501,6 +510,7 @@ def core(self, data):
         used_g = gOrig
         stroke = 0
         fine_tuned = False
+        coarse_limit_reached = False
         lng = False
 
         # SYMMETRIC APPROACH
@@ -529,6 +539,15 @@ def core(self, data):
                 v2 = 1.-abs((mult(Q.T, QOldF)).diagonal()).min(axis=0)
                 convergence_fine.append(v2)
 
+                if self.g != self.fine_g \
+                   and coarse_limit is not None \
+                   and convergence[round] < coarse_limit \
+                   and not coarse_limit_reached:
+                    if verbose:
+                        print 'Coarse convergence, switching to fine cost...'
+                    used_g = gFine
+                    coarse_limit_reached = True
+
                 if convergence[round] < limit:
                     if fine_tuning and (not fine_tuned):
                         if verbose:
@@ -569,7 +588,7 @@ def core(self, data):
                 # Show the progress...
                 if verbose:
                     msg = ('Step no. %d,'
-                           ' convergence: %.3f' % (round+1,convergence[round]))
+                           ' convergence: %.7f' % (round+1,convergence[round]))
                     print msg
 



On Thu, 15 Nov 2012, Kundu, Prantik (NIH/NIMH) [F] wrote:

> Hello All

> As you have guessed, the two-stage ica_nodes.py in AFNI is a custom modification, intended to make possible one contrast to bootstrap another, for example enabling the use of tanh contrast in 'somewhat' noisy data by bootstrapping convergence with the more robust pow3 contrast. If primary and secondary are both given as the same contrast function, then the secondary step is skipped and the standard MDP algo continues.  This modification is not in the MATLAB reference, and is something suited to the peculiarities of fMRI data. We know it improves sensitivity to spontaneous BOLD fluctuations in the course of denoising based on our NMR signal decay measures of BOLD signal quality, but otherwise it would be hard to tell the difference. Perhaps Tiziano would care to run this implementation on some reference data before deciding to integrate the patch into the mainline MDP. Importantly, the MDP as distributed in AFNI is isolated to the runtime directory of the code its distributed with (meica.libs), and is not made globally accessible.

> -Prantik
> ________________________________________
> From: Tiziano Zito [tiziano.zito at bccn-berlin.de]
> Sent: Thursday, November 15, 2012 3:52 PM
> To: NeuroDebian Development; Kundu, Prantik (NIH/NIMH) [F]
> Cc: mdp-dev
> Subject: Re: MDP patches in AFNI

> Hi,

> I had a quick look at the AFNI localized copy of MDP. It seems they
> imported a fairly recent git snapshot. If the ica_nodes.py
> patch would get integrated in MDP, I think they could use the
> released MDP 3.3 version without any incompatibility.

> Is this 2 stage convergence switching also present in the original
> MATLAB reference version of FastICA? I am far away from the ICA
> world since a long time now, so I may have missed some important
> developments lately. If this approach is indeed novel and introduced
> by PK, wouldn't be better to allow the user to choose which
> approach to use -- standard or 2-stage -- at instantiation time?
> This would make MDP users happy, because their results are not going
> to look different after upgrading to the next MDP release, and AFNI
> could still use exactly the same algorithm their are using now...

> What do you think?

> Ciao,
> Tiziano

> PS: Ccing the mdp-dev mailing list.

> On Thu 15 Nov, 15:08, Yaroslav Halchenko wrote:
> > Hi Tiziano,

> > NB Ccing PK -- original author of the patch

> > As you might (or not) know, AFNI now carries a copy of MDP (of state
> > somewhere before actual 3.3 release as far as I see) with slight
> > patching of ica_nodes.py to have 2 stages of convergence with switching
> > to gFine later in the loop (PK, please correct me if I am wrong)

> > I am attaching patch for your consideration to get absorbed into MDP
> > since we would like to avoid maintaining two copies of it ;)

> > if needed (may be there is more patches, PK?) you can find
> > complete AFNI source with this mdp copy at

> > http://git.debian.org/?p=pkg-exppsy/afni.git

> > branch upstream.

> > --
> > Yaroslav O. Halchenko
> > Postdoctoral Fellow,   Department of Psychological and Brain Sciences
> > Dartmouth College, 419 Moore Hall, Hinman Box 6207, Hanover, NH 03755
> > Phone: +1 (603) 646-9834                       Fax: +1 (603) 646-1419
> > WWW:   http://www.linkedin.com/in/yarik

> > --- pkundu/meica.libs/mdp/nodes/ica_nodes_old.py      2012-11-15 14:33:53.000000000 -0500
> > +++ pkundu/meica.libs/mdp/nodes/ica_nodes.py  2012-11-15 14:33:53.000000000 -0500
> > @@ -1,3 +1,4 @@
> > +
> >  __docformat__ = "restructuredtext en"

> >  import math
> > @@ -325,8 +326,8 @@
> >      def __init__(self, approach = 'defl', g = 'pow3', guess = None,
> >                   fine_g = 'pow3', mu = 1,
> >                   sample_size = 1, fine_tanh = 1, fine_gaus = 1,
> > -                 max_it = 1000, max_it_fine = 100,
> > -                 failures = 5, limit = 0.001, verbose = False,
> > +                 max_it = 5000, max_it_fine = 100,
> > +                 failures = 5, primary_limit=0.01, limit = 0.001,  verbose = False,
> >                   whitened = False, white_comp = None, white_parm = None,
> >                   input_dim = None, dtype=None):
> >          """
> > @@ -346,7 +347,9 @@
> >                        It is passed directly to the WhiteningNode constructor.
> >                        Ex: white_parm = { 'svd' : True }

> > -        limit -- convergence threshold.
> > +        limit -- final convergence threshold.
> > +
> > +     primary_limit -- initial convergence threshold, to switch to fine function, (i.e. linear to non-linear). PK 26-6-12.

> >          Specific for FastICA:

> > @@ -417,6 +420,7 @@
> >          self.fine_gaus = fine_gaus
> >          self.max_it = max_it
> >          self.max_it_fine = max_it_fine
> > +        self.primary_limit = primary_limit
> >          self.failures = failures
> >          self.guess = guess

> > @@ -458,6 +462,7 @@
> >                  guess = mult(guess, self.white.get_recmatrix(transposed=1))

> >          limit = self.limit
> > +        primary_limit = self.primary_limit
> >          max_it = self.max_it
> >          max_it_fine = self.max_it_fine
> >          failures = self.failures
> > @@ -501,6 +506,7 @@
> >          used_g = gOrig
> >          stroke = 0
> >          fine_tuned = False
> > +        in_secondary = False
> >          lng = False

> >          # SYMMETRIC APPROACH
> > @@ -529,10 +535,16 @@
> >                  v2 = 1.-abs((mult(Q.T, QOldF)).diagonal()).min(axis=0)
> >                  convergence_fine.append(v2)

> > +                if self.g!=self.fine_g and convergence[round] < primary_limit and not in_secondary:
> > +                    if verbose:
> > +                        print 'Primary convergence, switching to fine cost...'
> > +                    used_g = gFine
> > +                    in_secondary = True
> > +
> >                  if convergence[round] < limit:
> >                      if fine_tuning and (not fine_tuned):
> >                          if verbose:
> > -                            print 'Initial convergence, fine-tuning...'
> > +                            print 'Secondary convergence, fine-tuning...'
> >                          fine_tuned = True
> >                          used_g = gFine
> >                          mu = muK * self.mu
> > @@ -569,7 +581,7 @@
> >                  # Show the progress...
> >                  if verbose:
> >                      msg = ('Step no. %d,'
> > -                           ' convergence: %.3f' % (round+1,convergence[round]))
> > +                           ' convergence: %.7f' % (round+1,convergence[round]))
> >                      print msg




> _______________________________________________
> Neurodebian-devel mailing list
> Neurodebian-devel at lists.alioth.debian.org
> http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/neurodebian-devel


-- 
Yaroslav O. Halchenko
Postdoctoral Fellow,   Department of Psychological and Brain Sciences
Dartmouth College, 419 Moore Hall, Hinman Box 6207, Hanover, NH 03755
Phone: +1 (603) 646-9834                       Fax: +1 (603) 646-1419
WWW:   http://www.linkedin.com/in/yarik        



More information about the Neurodebian-devel mailing list