В этом примере показано, как использовать rica
распутывать смешанные звуковые сигналы. Можно использовать rica
выполнять независимый анализ компонента (ICA), когда предварительное отбеливание включено как шаг предварительной обработки. Модель ICA
Здесь,-by-1 вектор из смешанных сигналов,-by-1 вектор из значений смещения, смесительная матрица и-by-1 вектор из исходных сигналов. Предположим сначала, что это - квадратная матрица. Если вы знаете и, можно восстановить исходный сигнал с данных:
Используя rica
функция, можно выполнить это восстановление даже, не зная смесительной матрицы или среднего значения. Учитывая набор нескольких наблюдений..., rica
извлекает исходные сигналы....
Загрузите набор шести звуковых файлов, которые поставляются с MATLAB®. Обрежьте каждый файл к 10 000 выборок.
files = {'chirp.mat' 'gong.mat' 'handel.mat' 'laughter.mat' 'splat.mat' 'train.mat'}; S = zeros(10000,6); for i = 1:6 test = load(files{i}); y = test.y(1:10000,1); S(:,i) = y; end
Смешайте сигналы вместе при помощи случайной матрицы смешивания и добавьте случайное смещение.
rng default % For reproducibility mixdata = S*randn(6) + randn(1,6);
Чтобы слушать исходные звуки, выполните этот код:
for i = 1:6 disp(i); sound(S(:,i)); pause; end
Чтобы слушать смешанные звуки, выполните этот код:
for i = 1:6 disp(i); sound(mixdata(:,i)); pause; end
Постройте сигналы.
figure for i = 1:6 subplot(2,6,i) plot(S(:,i)) title(['Sound ',num2str(i)]) subplot(2,6,i+6) plot(mixdata(:,i)) title(['Mix ',num2str(i)]) end
Исходные сигналы имеют четкую структуру. Смешанные сигналы имеют намного меньше структуры.
Чтобы разделить сигналы эффективно, "предварительно побелите" сигналы при помощи prewhiten
функция, которая появляется в конце этого примера. Эта функция преобразовывает mixdata
так, чтобы это имело нулевое среднее значение и единичную ковариацию.
Идея следующая. Если нулевой средний источник со статистически независимыми компонентами, то
Затем среднее значение и ковариация
Предположим, что вы знаете и. На практике вы оценили бы эти количества от демонстрационного среднего значения и ковариации столбцов. Можно решить для в терминах
Последнее уравнение содержит, даже когда не квадратная обратимая матрица.
Предположим, что это - матрица левых собственных векторов положительной полуопределенной матрицы и является матрицей собственных значений. То
То
Существует много смесительных матриц, которые удовлетворяют этому последнему уравнению. Если ортонормированной матрицей, то
Замена в уравнение для,
предварительно побеленные данные. rica
вычисляет неизвестную матрицу под предположением, что компоненты максимально независимы.
mixdata = prewhiten(mixdata);
Супергауссов источник имеет резкий пик около нуля, такого как гистограмма звукового 1 показывает.
figure histogram(S(:,1))
Выполните ICA Реконструкции при выяснении шесть функций. Укажите, что каждый источник является супергауссовым.
q = 6;
Mdl = rica(mixdata,q,'NonGaussianityIndicator',ones(6,1));
Извлеките функции. Если несмесительная процедура успешна, функции пропорциональны исходным сигналам.
unmixed = transform(Mdl,mixdata);
Постройте исходные и несмешанные сигналы.
figure for i = 1:6 subplot(2,6,i) plot(S(:,i)) title(['Sound ',num2str(i)]) subplot(2,6,i+6) plot(unmixed(:,i)) title(['Unmix ',num2str(i)]) end
Порядок несмешанных сигналов отличается, чем первоначальный заказ. Переупорядочьте столбцы так, чтобы несмешанные сигналы совпадали с соответствующими исходными сигналами. Масштабируйте несмешанные сигналы иметь те же нормы как соответствующие исходные сигналы. (rica
не может идентифицировать шкалу исходных сигналов, потому что любая шкала может привести к той же смеси сигнала.)
unmixed = unmixed(:,[2,5,4,6,3,1]); for i = 1:6 unmixed(:,i) = unmixed(:,i)/norm(unmixed(:,i))*norm(S(:,i)); end
Постройте исходные и несмешанные сигналы.
figure for i = 1:6 subplot(2,6,i) plot(S(:,i)) ylim([-1,1]) title(['Sound ',num2str(i)]) subplot(2,6,i+6) plot(unmixed(:,i)) ylim([-1,1]) title(['Unmix ',num2str(i)]) end
Несмешанные сигналы выглядят похожими на исходные сигналы. Чтобы слушать несмешанные звуки, выполните этот код.
for i = 1:6 disp(i); sound(unmixed(:,i)); pause; end
Вот код для prewhiten
функция.
function Z = prewhiten(X) % X = N-by-P matrix for N observations and P predictors % Z = N-by-P prewhitened matrix % 1. Size of X. [N,P] = size(X); assert(N >= P); % 2. SVD of covariance of X. We could also use svd(X) to proceed but N % can be large and so we sacrifice some accuracy for speed. [U,Sig] = svd(cov(X)); Sig = diag(Sig); Sig = Sig(:)'; % 3. Figure out which values of Sig are non-zero. tol = eps(class(X)); idx = (Sig > max(Sig)*tol); assert(~all(idx == 0)); % 4. Get the non-zero elements of Sig and corresponding columns of U. Sig = Sig(idx); U = U(:,idx); % 5. Compute prewhitened data. mu = mean(X,1); Z = bsxfun(@minus,X,mu); Z = bsxfun(@times,Z*U,1./sqrt(Sig)); end
ReconstructionICA
| rica
| sparsefilt
| SparseFiltering