Python FastICA.get_mixing_matrix方法代码示例

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
def test_fastica_nowhiten():
    m = [[0, 1], [1, 0]]
    ica = FastICA(whiten=False, random_state=0)
    ica.fit(m)
    ica.get_mixing_matrix()

    # test for issue #697
    with warnings.catch_warnings(record=True) as w:
        warnings.simplefilter("always")
        ica = FastICA(n_components=1, whiten=False, random_state=0)
        ica.fit(m)  # should raise warning
        assert_true(len(w) == 1)  # 1 warning should be raised

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
def independent_component(x, y):
    clf = FastICA(random_state=1)
    clf.fit(x.reshape(-1, 1), y)
    comp = clf.components_[0][0]
    mm = clf.get_mixing_matrix()[0][0]
    sources = clf.sources_.flatten()
    src_max = max(sources)
    src_min = min(sources)
    return [comp, mm, src_max, src_min]

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
def independent_component( (x, y) ):
    lenX = len(x)
    newX = np.array(x).reshape(lenX, 1)
    g = FastICA()
    g.fit(newX, y)
    ret = [0.0, 0.0, 0.0, 0.0]
    ret[0] = g.components_[0][0]
    ret[1] = g.get_mixing_matrix()[0][0]
    sources = g.sources_.flatten()
    ret[2] = max(sources)
    ret[3] = min(sources)
    return ret

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
    def RunICAScikit(q):
      totalTimer = Timer()

      # Load input dataset.
      data = np.genfromtxt(self.dataset, delimiter=',')

      s = re.search('-s (\d+)', options)
      s = 0 if not s else int(s.group(1))

      try:
        # Perform ICA.
        with totalTimer:
          model = FastICA(random_state=s)
          ic = model.fit(data).transform(data)
          mixing = model.get_mixing_matrix()
      except Exception as e:
        q.put(-1)
        return -1

      time = totalTimer.ElapsedTime()
      q.put(time)
      return time

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]

#.........这里部分代码省略.........

    # functions of data
    minfnK    = lambda data: -kurt(data)**2
    minfnNEnt = lambda data: -negentropy(data)
    minfnLC10 = lambda data: -logcosh10(data)
    minfnLC15 = lambda data: -logcosh15(data)
    minfnLC20 = lambda data: -logcosh20(data)
    minfnNExp = lambda data: -negexp(data)

    # functions of the rotation angle, given W as the data
    minAngleFnK    = lambda theta: minfnK(dot(rotMat(theta)[0,:], W.T).T)
    minAngleFnNEnt = lambda theta: minfnNEnt(dot(rotMat(theta)[0,:], W.T).T)
    minAngleFnLC10 = lambda theta: minfnLC10(dot(rotMat(theta)[0,:], W.T).T)
    minAngleFnLC15 = lambda theta: minfnLC15(dot(rotMat(theta)[0,:], W.T).T)
    minAngleFnLC20 = lambda theta: minfnLC20(dot(rotMat(theta)[0,:], W.T).T)
    minAngleFnNExp = lambda theta: minfnNExp(dot(rotMat(theta)[0,:], W.T).T)

    #########
    # Chosen objective function. Change this line to change which objective is used.
    #########
    minDataFn = minfnK 

    minAngleFn = lambda theta: minDataFn(dot(rotMat(theta)[0,:], W.T).T)

    if doFastICA:
        # Use FastICA from sklearn
        #pdb.set_trace()
        from sklearn.decomposition import FastICA
        rng = random.RandomState(1)
        ica = FastICA(random_state = rng, whiten = False)
        ica.fit(W)
        Recon = ica.transform(W)  # Estimate the sources
        #S_fica /= S_fica.std(axis=0)   # (should already be done)
        Ropt = ica.get_mixing_matrix()
    else:
        # Manually fit angle using fmin_bfgs
        angle0 = 0
        xopt = fmin_bfgs(minAngleFn, angle0)
        xopt = xopt[0] % pi
        Ropt = rotMat(xopt)
        Recon = dot(W, Ropt.T)

    mnval = array([minAngleFn(aa) for aa in thetas])

    pyplot.figure(4)
    pyplot.title('objective vs. angle')
    #pyplot.plot(thetas, kurt1, 'bo-', thetas, mnval, 'k', xopt, minAngleFn(xopt), 'ko')
    pyplot.plot(thetas, mnval, 'b')
    if not doFastICA:
        pyplot.hold(True)
        pyplot.plot(xopt, minAngleFn(xopt), 'ko')

    pyplot.figure(5)
    pyplot.title('different gaussianness measures vs. angle')
    pyplot.subplot(6,1,1); pyplot.title('Kurt'); pyplot.plot(thetas, array([minAngleFnK(aa) for aa in thetas]))
    pyplot.subplot(6,1,2); pyplot.title('NegEnt'); pyplot.plot(thetas, array([minAngleFnNEnt(aa) for aa in thetas]))
    pyplot.subplot(6,1,3); pyplot.title('LogCosh10'); pyplot.plot(thetas, array([minAngleFnLC10(aa) for aa in thetas]))
    pyplot.subplot(6,1,4); pyplot.title('LogCosh15'); pyplot.plot(thetas, array([minAngleFnLC15(aa) for aa in thetas]))
    pyplot.subplot(6,1,5); pyplot.title('LogCosh20'); pyplot.plot(thetas, array([minAngleFnLC20(aa) for aa in thetas]))
    pyplot.subplot(6,1,6); pyplot.title('NegExp'); pyplot.plot(thetas, array([minAngleFnNExp(aa) for aa in thetas]))
    
    print 'kurt(r1) =', kurt(Recon[:,0])
    print 'kurt(r2) =', kurt(Recon[:,1])

    print
    print 'objective(s1) =', minDataFn(s1)

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
            if not x==y:  # skip diagonal
                plt.text(x-0.3,y,'%.2f'%cc_pcorr[x,y],size=6,color='red')

    plt.savefig(basedir+'9_correlation_analysis/pcorr_corrcoefs.pdf',format='pdf')


#do clustering

if 1==1:
    dst=pdist(data_pcorr[2:,:])
    Z=linkage(dst,method='complete')
    plt.figure(figsize=(14,12))
    dendrogram(Z,labels=tasknames_pcorr)
    plt.savefig(basedir+'9_correlation_analysis/pcorr_task_cluster.pdf',format='pdf')

# decompose connections using ICA and save adjacency matrices

data_pcorr_fmri=data_pcorr[2:,:]

if 1==0:
    ica = FastICA(n_components=20)
    S_ = ica.fit(data_pcorr_fmri.T).transform(data_pcorr_fmri.T)  # Get the estimated sources
    A_ = ica.get_mixing_matrix()  # Get estimated mixing matrix


                                  #ncomps=20
                                  #nmf=decomposition.ProjectedGradientNMF(n_components=ncomps,sparseness='components',init='nndsvd')
                                  #nmf.fit(data_pcorr_fmri+100)
#comps=nmf.components_

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
            # pl.quiver(x_axis, y_axis, x_axis, y_axis, zorder=11, width=0.01,
            pl.quiver(0, 0, x_axis, y_axis, zorder=11, width=0.01, scale=6,
                      color=color)

    pl.hlines(0, -3, 3)
    pl.vlines(0, -3, 3)
    pl.xlim(-3, 3)
    pl.ylim(-3, 3)
    pl.xlabel('x')
    pl.ylabel('y')

pl.subplot(2, 2, 1)
plot_samples(S / S.std())
pl.title('True Independent Sources')

axis_list = [pca.components_.T, ica.get_mixing_matrix()]
pl.subplot(2, 2, 2)
plot_samples(X / np.std(X), axis_list=axis_list)
pl.legend(['PCA', 'ICA'], loc='upper left')
pl.title('Observations')

pl.subplot(2, 2, 3)
plot_samples(S_pca_ / np.std(S_pca_, axis=0))
pl.title('PCA scores')

pl.subplot(2, 2, 4)
plot_samples(S_ica_ / np.std(S_ica_))
pl.title('ICA estimated sources')

pl.subplots_adjust(0.09, 0.04, 0.94, 0.94, 0.26, 0.26)

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
import numpy as np
import os
from sklearn.decomposition import FastICA

# Load the data.
data = np.genfromtxt('../data/radical.csv', delimiter=',')

# Peform FastICA
ica = FastICA()
S = ica.fit(data).transform(data)
A = ica.get_mixing_matrix() 

# show the results.
print S
print A

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
class ICA(object):
    """
    Wrapper for sklearn package.  Performs fast ICA (Independent Component Analysis)

    ICA has 4 methods:
       - fit(waveforms)
       update class instance with ICA fit

       - fit_transform()
       do what fit() does, but additionally return the projection onto ICA space

       - inverse_transform(A)
       inverses the decomposition, returns waveforms for an input A, using Z

       - get_params()
       returns metadata used for fits.
    """
    def __init__(self, num_components=10,
                 catalog_name='unknown',
                 whiten=True,
                 fun = 'logcosh',
                 fun_args = None,
                 max_iter = 600,
                 tol = .00001,
                 w_init = None,
                 random_state = None,
                 algorithm = 'parallel'):

        self._decomposition  = 'Fast ICA'
        self._num_components = num_components
        self._catalog_name   = catalog_name
        self._whiten         = whiten
        self._fun            = fun
        self._fun_args       = fun_args
        self._max_iter       = max_iter
        self._tol            = tol
        self._w_init         = w_init
        self._random_state   = random_state
        self._algorithm      = algorithm

        self._ICA = FastICA(n_components=self._num_components,
                             whiten       = self._whiten,
                             fun          = self._fun,
                             fun_args     = self._fun_args,
                             max_iter     = self._max_iter,
                             tol          = self._tol,
                             w_init       = self._w_init,
                             random_state = self._random_state,
                             algorithm    = self._algorithm)


    def fit(self,waveforms):
        # TODO make sure there are more columns than rows (transpose if not)
        # normalize waveforms
        self._waveforms = waveforms
        self._ICA.fit(self._waveforms)

    def fit_transform(self,waveforms):
        # TODO make sure there are more columns than rows (transpose if not)
        # normalize waveforms
        self._waveforms = waveforms
        self._A = self._ICA.fit_transform(self._waveforms)
        return self._A

    def inverse_transform(self,A):
        # convert basis back to waveforms using fit
        new_waveforms = self._ICA.inverse_transform(A)
        return new_waveforms

    def get_params(self):
        # TODO know what catalog was used! (include waveform metadata)
        params = self._ICA.get_params()
        params['num_components'] = params.pop('n_components')
        params['Decompositon'] = self._decomposition
        return params

    def get_basis(self):
        """ Return the ICA basis vectors (Z^\dagger)"""
        return self._ICA.get_mixing_matrix()

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
def test_fastica(add_noise=False):
    """ Test the FastICA algorithm on very simple data.
    """
    # scipy.stats uses the global RNG:
    rng = np.random.RandomState(0)
    n_samples = 1000
    # Generate two sources:
    s1 = (2 * np.sin(np.linspace(0, 100, n_samples)) > 0) - 1
    s2 = stats.t.rvs(1, size=n_samples)
    s = np.c_[s1, s2].T
    center_and_norm(s)
    s1, s2 = s

    # Mixing angle
    phi = 0.6
    mixing = np.array([[np.cos(phi),  np.sin(phi)],
                       [np.sin(phi), -np.cos(phi)]])
    m = np.dot(mixing, s)

    if add_noise:
        m += 0.1 * rng.randn(2, 1000)

    center_and_norm(m)

    algos = ['parallel', 'deflation']
    nls = ['logcosh', 'exp', 'cube']
    whitening = [True, False]
    for algo, nl, whiten in itertools.product(algos, nls, whitening):
        if whiten:
            k_, mixing_, s_ = fastica(m.T, fun=nl, algorithm=algo)
        else:
            X = PCA(n_components=2, whiten=True).fit_transform(m.T)
            k_, mixing_, s_ = fastica(X, fun=nl, algorithm=algo,
                                     whiten=False)
        s_ = s_.T
        # Check that the mixing model described in the docstring holds:
        if whiten:
            assert_almost_equal(s_, np.dot(np.dot(mixing_, k_), m))

        center_and_norm(s_)
        s1_, s2_ = s_
        # Check to see if the sources have been estimated
        # in the wrong order
        if abs(np.dot(s1_, s2)) > abs(np.dot(s1_, s1)):
            s2_, s1_ = s_
        s1_ *= np.sign(np.dot(s1_, s1))
        s2_ *= np.sign(np.dot(s2_, s2))

        # Check that we have estimated the original sources
        if add_noise == False:
            assert_almost_equal(np.dot(s1_, s1) / n_samples, 1, decimal=2)
            assert_almost_equal(np.dot(s2_, s2) / n_samples, 1, decimal=2)
        else:
            assert_almost_equal(np.dot(s1_, s1) / n_samples, 1, decimal=1)
            assert_almost_equal(np.dot(s2_, s2) / n_samples, 1, decimal=1)

    # Test FastICA class
    ica = FastICA(fun=nl, algorithm=algo, random_state=0)
    ica.fit(m.T)
    ica.get_mixing_matrix()
    assert_true(ica.components_.shape == (2, 2))
    assert_true(ica.sources_.shape == (1000, 2))

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]

#.........这里部分代码省略.........
        verbose : bool, str, int, or None
            If not None, override default verbose level (see mne.verbose).
            Defaults to self.verbose.

        Returns
        -------
        self : instance of ICA
            Returns the modified instance.
        """
        if self.current_fit != 'unfitted':
            raise RuntimeError('ICA decomposition has already been fitted. '
                               'Please start a new ICA session.')

        logger.info('Computing signal decomposition on raw data. '
                    'Please be patient, this may take some time')

        if picks is None:  # just use good data channels
            picks = pick_types(raw.info, meg=True, eeg=True, eog=False,
                               ecg=False, misc=False, stim=False,
                               exclude=raw.info['bads'])

        if self.max_n_components is None:
            self.max_n_components = len(picks)
            logger.info('Inferring max_n_components from picks.')

        self.ch_names = [raw.ch_names[k] for k in picks]

        data, self._pre_whitener = self._pre_whiten(raw[picks, start:stop][0],
                                                   raw.info, picks)

        to_ica, self._pca = self._prepare_pca(data, self.max_n_components)

        self._ica.fit(to_ica)
        self._mixing = self._ica.get_mixing_matrix().T
        self.current_fit = 'raw'

        return self

    @verbose
    def decompose_epochs(self, epochs, picks=None, verbose=None):
        """Run the ICA decomposition on epochs

        Parameters
        ----------
        epochs : instance of Epochs
            The epochs. The ICA is estimated on the concatenated epochs.
        picks : array-like
            Channels to be included relative to the channels already picked on
            epochs-initialization. This selection remains throughout the
            initialized ICA session.
        verbose : bool, str, int, or None
            If not None, override default verbose level (see mne.verbose).
            Defaults to self.verbose.

        Returns
        -------
        self : instance of ICA
            Returns the modified instance.
        """
        if self.current_fit != 'unfitted':
            raise RuntimeError('ICA decomposition has already been fitted. '
                               'Please start a new ICA session.')

        logger.info('Computing signal decomposition on epochs. '
                    'Please be patient, this may take some time')

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
def clusterize_dirichlet(*args, **kwargs):
    # TODO plotting when the classifier is NOT learned in the PCA(2-best) space (ellipses are wrong)
    """ Clustering and plotting with Dirichlet process GMM """
    ### Clustering
    try:
        from sklearn import mixture
        from scipy import linalg
        import pylab as pl
        import matplotlib as mpl
        from sklearn.decomposition import PCA, FastICA
    except:
        print "You need SciPy and scikit-learn"
        sys.exit(-1)

    models = []
    for arg in args:
        if apply_ica:
            for featurenb in range(len(arg[0])):
                if sum(arg[:, featurenb]) == 0.0:
                    arg = arg[:, range(featurenb+1)+range(featurenb+1, len(arg[0]))]
        if kwargs.get('em_gmm', False):
            dpgmm = mixture.GMM(n_components = 4, cvtype='full')
        else:
            dpgmm = mixture.DPGMM(n_components = 100, cvtype='full', alpha=1000.0)
        if kwargs.get('clf_on_pca', False):
            pca = PCA(2)
            dpgmm.fit(pca.fit(arg).transform(arg))
        else:
            dpgmm.fit(arg)
        print dpgmm
        models.append(copy.deepcopy(dpgmm))
        print raw_input("press any key to pass")

    ### Plotting
    color_iter = itertools.cycle (['r', 'g', 'b', 'c', 'm'])
    for i, (clf, data) in enumerate(zip(models, args)):
        if apply_ica:
            ica = FastICA(2)
            X_r = ica.fit(data).transform(data)
            print ica.get_mixing_matrix()
        else:
            pca = PCA(2)
            X_r = pca.fit(data).transform(data)
        print data
        print X_r
        print raw_input("press any key to pass")
        splot = pl.subplot((len(args)+1)/ 2, 2, 1+i)
        pl.scatter(X_r[:,0], X_r[:,1])
        if kwargs.get('clf_on_pca', False):
            Y_ = clf.predict(X_r)
        else:
            Y_ = clf.predict(data)
        for i, (mean, covar, color) in enumerate(zip(clf.means, clf.covars,
                                                     color_iter)):
            v, w = linalg.eigh(covar)
            u = w[0] / linalg.norm(w[0])
            # as the DP will not use every component it has access to
            # unless it needs it, we shouldn't plot the redundant
            # components.
            if not np.any(Y_ == i):
                continue
            #pl.scatter(data[Y_== i, 0], data[Y_== i, 1], .8, color=color)
            pl.scatter(X_r[Y_== i, 0], X_r[Y_== i, 1], .8, color=color)
            # Plot an ellipse to show the Gaussian component
            angle = np.arctan(u[1]/u[0])
            angle = 180 * angle / np.pi # convert to degrees
            ell = mpl.patches.Ellipse(mean, v[0], v[1], 180 + angle, color=color)
            ell.set_clip_box(splot.bbox)
            ell.set_alpha(0.5)
            splot.add_artist(ell)
        pl.xlim(X_r[:, 0].min(), X_r[:, 0].max())
        pl.ylim(X_r[:, 1].min(), X_r[:, 1].max())
        pl.xticks(())
        pl.yticks(())
        pl.title("Dirichlet process GMM")
    pl.show()

# 需要导入模块: from sklearn.decomposition import FastICA [as 别名]
# 或者: from sklearn.decomposition.FastICA import get_mixing_matrix [as 别名]
def test_fastica_simple(add_noise=False):
    """ Test the FastICA algorithm on very simple data.
    """
    rng = np.random.RandomState(0)
    # scipy.stats uses the global RNG:
    np.random.seed(0)
    n_samples = 1000
    # Generate two sources:
    s1 = (2 * np.sin(np.linspace(0, 100, n_samples)) > 0) - 1
    s2 = stats.t.rvs(1, size=n_samples)
    s = np.c_[s1, s2].T
    center_and_norm(s)
    s1, s2 = s

    # Mixing angle
    phi = 0.6
    mixing = np.array([[np.cos(phi), np.sin(phi)], [np.sin(phi), -np.cos(phi)]])
    m = np.dot(mixing, s)

    if add_noise:
        m += 0.1 * rng.randn(2, 1000)

    center_and_norm(m)

    # function as fun arg
    def g_test(x):
        return x ** 3, 3 * x ** 2

    algos = ["parallel", "deflation"]
    nls = ["logcosh", "exp", "cube", g_test]
    whitening = [True, False]
    for algo, nl, whiten in itertools.product(algos, nls, whitening):
        if whiten:
            k_, mixing_, s_ = fastica(m.T, fun=nl, algorithm=algo)
            assert_raises(ValueError, fastica, m.T, fun=np.tanh, algorithm=algo)
        else:
            X = PCA(n_components=2, whiten=True).fit_transform(m.T)
            k_, mixing_, s_ = fastica(X, fun=nl, algorithm=algo, whiten=False)
            assert_raises(ValueError, fastica, X, fun=np.tanh, algorithm=algo)
        s_ = s_.T
        # Check that the mixing model described in the docstring holds:
        if whiten:
            assert_almost_equal(s_, np.dot(np.dot(mixing_, k_), m))

        center_and_norm(s_)
        s1_, s2_ = s_
        # Check to see if the sources have been estimated
        # in the wrong order
        if abs(np.dot(s1_, s2)) > abs(np.dot(s1_, s1)):
            s2_, s1_ = s_
        s1_ *= np.sign(np.dot(s1_, s1))
        s2_ *= np.sign(np.dot(s2_, s2))

        # Check that we have estimated the original sources
        if not add_noise:
            assert_almost_equal(np.dot(s1_, s1) / n_samples, 1, decimal=2)
            assert_almost_equal(np.dot(s2_, s2) / n_samples, 1, decimal=2)
        else:
            assert_almost_equal(np.dot(s1_, s1) / n_samples, 1, decimal=1)
            assert_almost_equal(np.dot(s2_, s2) / n_samples, 1, decimal=1)

    # Test FastICA class
    ica = FastICA(fun=nl, algorithm=algo, random_state=0)
    ica.fit(m.T)
    ica.get_mixing_matrix()
    assert_true(ica.components_.shape == (2, 2))
    assert_true(ica.sources_.shape == (1000, 2))

    for fn in [np.tanh, "exp(-.5(x^2))"]:
        ica = FastICA(fun=fn, algorithm=algo, random_state=0)
        assert_raises(ValueError, ica.fit, m.T)

    assert_raises(TypeError, FastICA(fun=xrange(10)).fit, m.T)