Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
added UWedgeSep and python modular examples for all current ICA methods
- Loading branch information
1 parent
25f81d5
commit 53014e3
Showing
10 changed files
with
525 additions
and
2 deletions.
There are no files selected for viewing
57 changes: 57 additions & 0 deletions
57
examples/undocumented/python_modular/graphical/converter_ffsep_bss.py
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,57 @@ | ||
""" | ||
Blind Source Separation using the FFSep Algorithm with Shogun | ||
Based on the example from scikit-learn | ||
http://scikit-learn.org/ | ||
Kevin Hughes 2013 | ||
""" | ||
|
||
|
||
import numpy as np | ||
import pylab as pl | ||
|
||
from shogun.Features import RealFeatures | ||
from shogun.Converter import FFSep | ||
|
||
# Generate sample data | ||
np.random.seed(0) | ||
n_samples = 2000 | ||
time = np.linspace(0, 10, n_samples) | ||
|
||
# Source Signals | ||
s1 = np.sin(2 * time) # sin wave | ||
s2 = np.sign(np.sin(3 * time)) # square wave | ||
S = np.c_[s1, s2] | ||
S += 0.2 * np.random.normal(size=S.shape) # add noise | ||
|
||
# Standardize data | ||
S /= S.std(axis=0) | ||
S = S.T | ||
|
||
# Mixing Matrix | ||
A = np.array([[1, 0.5], [0.5, 1]]) | ||
|
||
# Mix Signals | ||
X = np.dot(A,S) | ||
mixed_signals = RealFeatures(X) | ||
|
||
# Separating | ||
ffsep = FFSep() | ||
signals = ffsep.apply(mixed_signals) | ||
S_ = signals.get_feature_matrix() | ||
A_ = ffsep.get_mixing_matrix(); | ||
|
||
# Plot results | ||
pl.figure() | ||
pl.subplot(3, 1, 1) | ||
pl.plot(S.T) | ||
pl.title('True Sources') | ||
pl.subplot(3, 1, 2) | ||
pl.plot(X.T) | ||
pl.title('Mixed Sources') | ||
pl.subplot(3, 1, 3) | ||
pl.plot(S_.T) | ||
pl.title('Estimated Sources') | ||
pl.subplots_adjust(0.09, 0.04, 0.94, 0.94, 0.26, 0.36) | ||
pl.show() |
59 changes: 59 additions & 0 deletions
59
examples/undocumented/python_modular/graphical/converter_jedi_bss.py
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,59 @@ | ||
""" | ||
Blind Source Separation using the JediSep Algorithm with Shogun | ||
Based on the example from scikit-learn | ||
http://scikit-learn.org/ | ||
Kevin Hughes 2013 | ||
""" | ||
|
||
|
||
import numpy as np | ||
import pylab as pl | ||
|
||
from shogun.Features import RealFeatures | ||
from shogun.Converter import JediSep | ||
|
||
# Generate sample data | ||
np.random.seed(0) | ||
n_samples = 2000 | ||
time = np.linspace(0, 10, n_samples) | ||
|
||
# Source Signals | ||
s1 = np.sin(2 * time) # sin wave | ||
s2 = np.sign(np.sin(3 * time)) # square wave | ||
S = np.c_[s1, s2] | ||
S += 0.2 * np.random.normal(size=S.shape) # add noise | ||
|
||
# Standardize data | ||
S /= S.std(axis=0) | ||
S = S.T | ||
|
||
# Mixing Matrix | ||
A = np.array([[1, 0.5], [0.5, 1]]) | ||
|
||
# Mix Signals | ||
X = np.dot(A,S) | ||
mixed_signals = RealFeatures(X) | ||
|
||
# Separating | ||
jedi = JediSep() | ||
signals = jedi.apply(mixed_signals) | ||
S_ = signals.get_feature_matrix() | ||
A_ = jedi.get_mixing_matrix(); | ||
|
||
print A_ | ||
|
||
# Plot results | ||
pl.figure() | ||
pl.subplot(3, 1, 1) | ||
pl.plot(S.T) | ||
pl.title('True Sources') | ||
pl.subplot(3, 1, 2) | ||
pl.plot(X.T) | ||
pl.title('Mixed Sources') | ||
pl.subplot(3, 1, 3) | ||
pl.plot(S_.T) | ||
pl.title('Estimated Sources') | ||
pl.subplots_adjust(0.09, 0.04, 0.94, 0.94, 0.26, 0.36) | ||
pl.show() |
57 changes: 57 additions & 0 deletions
57
examples/undocumented/python_modular/graphical/converter_sobi_bss.py
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,57 @@ | ||
""" | ||
Blind Source Separation using the SOBI Algorithm with Shogun | ||
Based on the example from scikit-learn | ||
http://scikit-learn.org/ | ||
Kevin Hughes 2013 | ||
""" | ||
|
||
|
||
import numpy as np | ||
import pylab as pl | ||
|
||
from shogun.Features import RealFeatures | ||
from shogun.Converter import SOBI | ||
|
||
# Generate sample data | ||
np.random.seed(0) | ||
n_samples = 2000 | ||
time = np.linspace(0, 10, n_samples) | ||
|
||
# Source Signals | ||
s1 = np.sin(2 * time) # sin wave | ||
s2 = np.sign(np.sin(3 * time)) # square wave | ||
S = np.c_[s1, s2] | ||
S += 0.2 * np.random.normal(size=S.shape) # add noise | ||
|
||
# Standardize data | ||
S /= S.std(axis=0) | ||
S = S.T | ||
|
||
# Mixing Matrix | ||
A = np.array([[1, 0.5], [0.5, 1]]) | ||
|
||
# Mix Signals | ||
X = np.dot(A,S) | ||
mixed_signals = RealFeatures(X) | ||
|
||
# Separating | ||
sobi = SOBI() | ||
signals = sobi.apply(mixed_signals) | ||
S_ = signals.get_feature_matrix() | ||
A_ = sobi.get_mixing_matrix(); | ||
|
||
# Plot results | ||
pl.figure() | ||
pl.subplot(3, 1, 1) | ||
pl.plot(S.T) | ||
pl.title('True Sources') | ||
pl.subplot(3, 1, 2) | ||
pl.plot(X.T) | ||
pl.title('Mixed Sources') | ||
pl.subplot(3, 1, 3) | ||
pl.plot(S_.T) | ||
pl.title('Estimated Sources') | ||
pl.subplots_adjust(0.09, 0.04, 0.94, 0.94, 0.26, 0.36) | ||
pl.show() |
56 changes: 56 additions & 0 deletions
56
examples/undocumented/python_modular/graphical/converter_uwedge_bss.py
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,56 @@ | ||
""" | ||
Blind Source Separation using the UWedgeSep Algorithm with Shogun | ||
Kevin Hughes 2013 | ||
""" | ||
|
||
|
||
import numpy as np | ||
import pylab as pl | ||
|
||
from shogun.Features import RealFeatures | ||
from shogun.Converter import UWedgeSep | ||
|
||
# Generate sample data | ||
np.random.seed(0) | ||
n_samples = 2000 | ||
time = np.linspace(0, 1, n_samples) | ||
|
||
# Source Signals | ||
s1 = np.sin(2*3.14*55*time) # sin wave | ||
s2 = np.cos(2*3.14*100*time) # cos wave | ||
S = np.c_[s1, s2] | ||
S += 0.1 * np.random.normal(size=S.shape) # add noise | ||
|
||
# Standardize data | ||
S /= S.std(axis=0) | ||
S = S.T | ||
|
||
# Mixing Matrix | ||
A = np.array([[1, 0.85], [0.55, 1]]) | ||
#print A | ||
|
||
# Mix Signals | ||
X = np.dot(A,S) | ||
mixed_signals = RealFeatures(X) | ||
|
||
# Separating | ||
uwedge = UWedgeSep() | ||
signals = uwedge.apply(mixed_signals) | ||
S_ = signals.get_feature_matrix() | ||
A_ = uwedge.get_mixing_matrix(); | ||
#print A_ | ||
|
||
# Plot results | ||
pl.figure() | ||
pl.subplot(3, 1, 1) | ||
pl.plot(S.T) | ||
pl.title('True Sources') | ||
pl.subplot(3, 1, 2) | ||
pl.plot(X.T) | ||
pl.title('Mixed Sources') | ||
pl.subplot(3, 1, 3) | ||
pl.plot(S_.T) | ||
pl.title('Estimated Sources') | ||
pl.subplots_adjust(0.09, 0.04, 0.94, 0.94, 0.26, 0.36) | ||
pl.show() |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,140 @@ | ||
/* | ||
* This program is free software; you can redistribute it and/or modify | ||
* it under the terms of the GNU General Public License as published by | ||
* the Free Software Foundation; either version 3 of the License, or | ||
* (at your option) any later version. | ||
* | ||
* Written (W) 2013 Kevin Hughes | ||
*/ | ||
|
||
#include <shogun/converter/ica/UWedgeSep.h> | ||
|
||
#include <shogun/features/DenseFeatures.h> | ||
|
||
#ifdef HAVE_EIGEN3 | ||
|
||
#include <shogun/mathematics/Math.h> | ||
#include <shogun/mathematics/eigen3.h> | ||
#include <shogun/mathematics/ajd/UWedge.h> | ||
|
||
using namespace Eigen; | ||
|
||
typedef Matrix< float64_t, Dynamic, 1, ColMajor > EVector; | ||
typedef Matrix< float64_t, Dynamic, Dynamic, ColMajor > EMatrix; | ||
|
||
using namespace shogun; | ||
|
||
namespace { EMatrix cor(EMatrix x, int tau = 0, bool mean_flag = true); }; | ||
|
||
CUWedgeSep::CUWedgeSep() : CConverter() | ||
{ | ||
m_tau = SGVector<float64_t>(4); | ||
m_tau[0]=0; m_tau[1]=1; m_tau[2]=2; m_tau[3]=3; | ||
|
||
init(); | ||
} | ||
|
||
void CUWedgeSep::init() | ||
{ | ||
m_mixing_matrix = SGMatrix<float64_t>(); | ||
|
||
SG_ADD(&m_tau, "tau", "tau vector", MS_AVAILABLE); | ||
SG_ADD(&m_mixing_matrix, "mixing_matrix", "m_mixing_matrix", MS_NOT_AVAILABLE); | ||
} | ||
|
||
CUWedgeSep::~CUWedgeSep() | ||
{ | ||
} | ||
|
||
void CUWedgeSep::set_tau(SGVector<float64_t> tau) | ||
{ | ||
m_tau = tau; | ||
} | ||
|
||
SGVector<float64_t> CUWedgeSep::get_tau() const | ||
{ | ||
return m_tau; | ||
} | ||
|
||
SGMatrix<float64_t> CUWedgeSep::get_mixing_matrix() const | ||
{ | ||
return m_mixing_matrix; | ||
} | ||
|
||
CFeatures* CUWedgeSep::apply(CFeatures* features) | ||
{ | ||
REQUIRE(features, "features is null"); | ||
SG_REF(features); | ||
|
||
SGMatrix<float64_t> X = ((CDenseFeatures<float64_t>*)features)->get_feature_matrix(); | ||
|
||
int n = X.num_rows; | ||
int m = X.num_cols; | ||
int N = m_tau.vlen; | ||
|
||
Eigen::Map<EMatrix> EX(X.matrix,n,m); | ||
|
||
// Compute Correlation Matrices | ||
index_t * M_dims = SG_MALLOC(index_t, 3); | ||
M_dims[0] = n; | ||
M_dims[1] = n; | ||
M_dims[2] = N; | ||
SGNDArray< float64_t > M(M_dims, 3); | ||
|
||
for (int t = 0; t < N; t++) | ||
{ | ||
Eigen::Map<EMatrix> EM(M.get_matrix(t),n,n); | ||
EM = cor(EX,m_tau[t]); | ||
} | ||
|
||
// Diagonalize | ||
SGMatrix<float64_t> Q = CUWedge::diagonalize(M); | ||
Eigen::Map<EMatrix> EQ(Q.matrix,n,n); | ||
|
||
// Compute Mixing Matrix | ||
m_mixing_matrix = SGMatrix<float64_t>(n,n); | ||
Eigen::Map<EMatrix> C(m_mixing_matrix.matrix,n,n); | ||
C = EQ.inverse(); | ||
|
||
// Normalize Estimated Mixing Matrix | ||
for (int t = 0; t < C.cols(); t++) | ||
C.col(t) /= C.col(t).maxCoeff(); | ||
|
||
// Unmix | ||
EX = C.inverse() * EX; | ||
|
||
return features; | ||
} | ||
|
||
// Computing time delayed correlation matrix | ||
namespace | ||
{ | ||
EMatrix cor(EMatrix x, int tau, bool mean_flag) | ||
{ | ||
int m = x.rows(); | ||
int n = x.cols(); | ||
|
||
// Center the data | ||
if ( mean_flag ) | ||
{ | ||
EVector mean = x.rowwise().sum(); | ||
mean /= n; | ||
x = x.colwise() - mean; | ||
} | ||
|
||
// Time-delayed Signal Matrix | ||
EMatrix L = x.leftCols(n-tau); | ||
EMatrix R = x.rightCols(n-tau); | ||
|
||
// Compute Correlations | ||
EMatrix K(m,m); | ||
K = (L * R.transpose()) / (n-tau); | ||
|
||
// Symmetrize | ||
K = (K + K.transpose()) / 2.0; | ||
|
||
return K; | ||
} | ||
}; | ||
|
||
#endif // HAVE_EIGEN3 |
Oops, something went wrong.