ZF vs MMSE
Пространственное мультиплексирование 2×2 MIMO: сравнение детекторов ZF и MMSE
Пример демонстрирует работу системы пространственного мультиплексирования (Spatial Multiplexing) в канале MIMO 2×2 с рэлеевскими замираниями. По двум антеннам передаются два независимых QPSK-потока, а на приёмной стороне применяется линейное разделение потоков. Сравниваются два классических линейных детектора — Zero-Forcing (ZF) и Minimum Mean Square Error (MMSE) — по критериям MSE и BER, с визуализацией созвездий.
Модель приёма имеет вид:
где — вектор передаваемых символов, — матрица канала (), — комплексный белый гауссовский шум.
Zero-Forcing (ZF) полностью подавляет межпоточные искажения, обращая канал:
Недостаток: при малом сингулярном числе шум усиливается на коэффициент , что резко ухудшает MSE.
MMSE минимизирует среднеквадратичную ошибку , учитывая дисперсию шума :
Регуляризатор предотвращает взрывной рост шума в направлении слабых собственных чисел ценой небольшого смещения оценки (bias). Именно поэтому при низких и умеренных SNR и плохой обусловленности канала MMSE существенно выигрывает у ZF.
Реализация
Скрипт состоит из следующих этапов:
- Генерация данных — формируется матрица символов QPSK (2 потока × N символов) и случайная матрица канала . Канал искусственно обусловливается (через умножение на фиксированные матрицы), чтобы наглядно проявить преимущество MMSE.
- Канал и шум — сигнал проходит через , после чего к нему добавляется AWGN с помощью системного объекта
EngeeComms.AWGN(параметрыEbNo,BitsPerSymbol,SignalPower). - Детектирование — вычисляются матрицы и через устойчивое LU-разложение (оператор
\), и применяются к принятому сигналу. - Оценка качества — рассчитываются MSE для обоих детекторов и выводится число обусловленности .
- Визуализация — строятся 6 созвездий: принятый сигнал на двух антеннах (смесь потоков), результаты ZF и MMSE для каждого потока.
- BER-свип — в цикле по SNR методом Монте-Карло оценивается BER, а подсчёт ошибок выполняется через системный объект
EngeeComms.ErrorRateCalculation.
Использование системных объектов Engee (EngeeComms.AWGN, EngeeComms.ErrorRateCalculation) обеспечивает корректный расчёт энергетических параметров и статистики ошибок в соответствии со стандартами библиотеки, избавляя от рутинных операций и делая код пригодным для переноса в полноценные модели связи.
# Pkg.add(["LinearAlgebra", "Random", "Statistics"])
using LinearAlgebra, Random, Statistics, EngeeComms
Для начала выполним инициализации системы:
Random.seed!(2024)— фиксация ГСЧ для воспроизводимости.- Параметры:
N=4000символов/поток,SNR=8 дБ, мощность сигнала 1. Пересчёт для блока AWGN и для MMSE. - QPSK: созвездие нормировано на единичную мощность, формируется матрица символов .
- Канал : рэлеевская матрица искусственно обусловливается умножением на фиксированные матрицы, чтобы — это наглядно проявит преимущество MMSE над ZF.
y_ideal = H * s— сигнал проходит через канал (без шума): потоки смешиваются, и задача детекторов — их разделить.
Random.seed!(2024)
N = 4000
SNR_dB = 8.0
bits_per_symbol = 2
SignalPower = 1.0
EbNo_dB = SNR_dB - 10 * log10(bits_per_symbol)
σ² = 1.0 / (10^(SNR_dB / 10))
qpsk = ComplexF64[1+1im, -1+1im, -1-1im, 1-1im] ./ sqrt(2)
s = qpsk[rand(1:4, 2, N)]
H = (randn(2,2) .+ 1im .* randn(2,2)) ./ sqrt(2)
H = H * [1.0 0.3; 0.0 0.15] * [1.0 0.0; 0.4 1.0]
y_ideal = H * s
Далее добавляем шум:
EngeeComms.AWGN(...)— создание системного объекта, моделирующего аддитивный белый гауссовский шум. Параметры:EbNo(отношение энергии бита к шуму),BitsPerSymbol=2(для QPSK),SignalPower=1.0(мощность сигнала).awgn_channel(vec(y_ideal))— прохождение сигнала через канал с шумом. Функцияvec()преобразует матрицу в вектор-столбец, так как API Engee принимает одномерные массивы.reshape(y_vec, size(y_ideal))— восстановление исходной формы матрицы для дальнейшей обработки детекторами.
awgn_channel = EngeeComms.AWGN(
EbNo = EbNo_dB,
BitsPerSymbol = bits_per_symbol,
SignalPower = SignalPower
)
y_vec = awgn_channel(vec(y_ideal))
y = reshape(y_vec, size(y_ideal))
Следующий этап — детектирование и оценка качества. Здесь происходит:
- Расчёт матриц детекторов
W_zfиW_mmseпо формулам псевдообращения и Винера-Хопфа. - Применение детекторов к принятому сигналу
y— разделение смешанных потоков. - Вычисление MSE — количественная оценка точности восстановления символов для обоих методов.
W_zf = (H' * H) \ H'
ŝ_zf = W_zf * y
W_mmse = (H' * H + σ² * I) \ H'
ŝ_mmse = W_mmse * y
mse_zf = mean(abs2, ŝ_zf - s)
mse_mmse = mean(abs2, ŝ_mmse - s)
println("MIMO 2×2: ZF vs MMSE")
println("SNR (Es/N0) = $SNR_dB дБ | Eb/No = $(round(EbNo_dB, digits=2)) дБ")
println("cond(H) = $(round(cond(H), digits=2))")
println("MSE ZF = $(round(mse_zf, digits=5))")
println("MSE MMSE = $(round(mse_mmse, digits=5))")
После этого идёт визуализация созвездий
draw(x; ...)— вспомогательная функция, строящая созвездие в два слоя:- полупрозрачное облако принятых точек (I/Q);
- красные кресты — идеальные позиции QPSK для сравнения.
- 6 графиков: принятый сигнал на двух антеннах (серые, смесь потоков), результаты ZF (синие) и MMSE (зелёные) для каждого из двух потоков.
function draw(x; title="", color=:blue)
p = scatter(real.(x), imag.(x);
markersize=3, alpha=0.35, color=color, title=title,
xlabel="I", ylabel="Q",
xlims=(-2.5, 2.5), ylims=(-2.5, 2.5),
aspect_ratio=:equal, legend=:topright, grid=true)
scatter!(p, real.(qpsk), imag.(qpsk);
markersize=8, color=:red, marker=:xcross)
return p
end
p_y1 = draw(y[1,:]; title="Rx антенна 1", color=:gray)
p_y2 = draw(y[2,:]; title="Rx антенна 2", color=:gray)
p_z1 = draw(ŝ_zf[1,:]; title="Поток 1 после ZF", color=:blue)
p_z2 = draw(ŝ_zf[2,:]; title="Поток 2 после ZF", color=:blue)
p_m1 = draw(ŝ_mmse[1,:]; title="Поток 1 после MMSE", color=:green)
p_m2 = draw(ŝ_mmse[2,:]; title="Поток 2 после MMSE", color=:green)
plot(p_y1, p_y2, p_z1, p_z2, p_m1, p_m2;
layout=(3,2), size=(900, 1100),
plot_title="MIMO 2×2: ZF vs MMSE (SNR = $SNR_dB дБ)",
left_margin=8Plots.mm, bottom_margin=8Plots.mm)
Вывод по графикам:
Верхний ряд — на антеннах видны 4 кластера (смесь потоков), на антенне 2 более размытые.
Средний ряд (ZF) — катастрофический разброс:
- Поток 1: огромное синее облако, точки заполняют всё поле
- Поток 2: ещё хуже — точки равномерно размазаны, структура QPSK полностью потеряна
Нижний ряд (MMSE) — кардинально лучше:
- Поток 1: 4 чётких зелёных кластера, идеальное восстановление
- Поток 2: один кластер слился в центре.
Инициализация ErrorRateCalculation
snrs_dB = 0:2:20— диапазон отношений сигнал/шум (от 0 до 20 дБ с шагом 2) для построения кривых BER.ber_zf, ber_mmse— пустые массивы для накопления результатов BER по каждому SNR для детекторов ZF и MMSE.Nsym = 20_000— количество символов в одном прогоне (итерации) для обеспечения статистической достоверности.n_iterations = 50— число итераций Монте-Карло (различных случайных реализаций канала) на каждую точку SNR.err_calc = EngeeComms.ErrorRateCalculation()— создание системного объекта Engee, который будет использоваться для подсчёта количества ошибок на последующих шагах цикла.
snrs_dB = 0:2:20
ber_zf, ber_mmse = Float64[], Float64[]
Nsym = 20_000
n_iterations = 50
err_calc = EngeeComms.ErrorRateCalculation()
В коде далее имеются два вложенных цикла. Внешний перебирает значения SNR от 0 до 20 дБ с шагом 2, инициализирует счётчики ошибок и создаёт объект EngeeComms.AWGN с пересчитанным для текущего шага значением . Внутренний цикл Монте-Карло выполняет 50 независимых реализаций: на каждой генерируется новая случайная матрица рэлеевского канала и блок передаваемых символов , сигнал проходит через канал с шумом, после чего применяются детекторы ZF и MMSE. Жёсткое решение относит каждый восстановленный символ к ближайшей точке QPSK, а системный объект ErrorRateCalculation подсчитывает количество ошибок, которые накапливаются по всем итерациям.
По завершении внутреннего цикла вычисляется SER как отношение суммарного числа ошибок к общему количеству переданных символов, а затем пересчитывается в BER с учётом того, что для QPSK с кодированием Грея . Полученные значения добавляются в массивы ber_zf и ber_mmse.
Финальный график строится в логарифмической шкале по оси ординат с двумя кривыми — ZF и MMSE.
for snr_db in snrs_dB
err_zf_total = 0
err_mmse_total = 0
ebno_db_sweep = snr_db - 10 * log10(bits_per_symbol)
awgn_sweep = EngeeComms.AWGN(EbNo = ebno_db_sweep, BitsPerSymbol = bits_per_symbol, SignalPower = SignalPower)
for _ in 1:n_iterations
H_ = (randn(2,2) .+ 1im.*randn(2,2)) ./ sqrt(2)
s_ = qpsk[rand(1:4, 2, Nsym)]
y_ideal_sweep = H_ * s_
y_sweep_vec = awgn_sweep(vec(y_ideal_sweep))
y_sweep = reshape(y_sweep_vec, size(y_ideal_sweep))
σ2_sweep = 1.0 / (10^(snr_db / 10))
ŝz = ((H_'*H_) \ H_') * y_sweep
ŝm = ((H_'*H_ + σ2_sweep*I) \ H_') * y_sweep
decide(x) = sign(real(x))/sqrt(2) + 1im*sign(imag(x))/sqrt(2)
bits_zf = decide.(ŝz)
bits_mmse = decide.(ŝm)
stats_zf = err_calc(vec(bits_zf), vec(s_))
stats_mmse = err_calc(vec(bits_mmse), vec(s_))
err_zf_total += stats_zf[2]
err_mmse_total += stats_mmse[2]
end
ser_zf = err_zf_total / (n_iterations * Nsym)
ser_mmse = err_mmse_total / (n_iterations * Nsym)
push!(ber_zf, ser_zf / 2.0)
push!(ber_mmse, ser_mmse / 2.0)
end
p_ber = plot(snrs_dB, ber_zf; yscale=:log10, marker=:circle, lw=2, label="ZF",
xlabel="SNR (Es/N0), дБ", ylabel="BER", title="BER: ZF vs MMSE (2×2 MIMO)",
grid=true, legend=:bottomleft)
plot!(p_ber, snrs_dB, ber_mmse; marker=:square, lw=2, label="MMSE", color=:green)
Почему графики совпали. В отличие от демонстрации созвездий, где канал был искусственно обусловлен (), в BER используется «чистый» Rayleigh-канал. Для случайной матрицы с комплексными гауссовскими элементами вероятность получить плохо обусловленную реализацию мала — в среднем невелико. При хорошей обусловленности и умеренных SNR регуляризатор в MMSE становится пренебрежимо мал, и оба детектора фактически сводятся к псевдообращению . Поэтому на графике кривые ZF и MMSE практически сливаются: теоретический выигрыш MMSE здесь «размазывается» по множеству благоприятных реализаций канала и не виден на фоне статистического усреднения. Чтобы выявить различие на BER-кривой, потребовалось бы либо зафиксировать плохо обусловленный канал (как в первой части примера), либо сдвинуть диапазон SNR в область очень низких значений.
Вывод
Пример наглядно демонстрирует фундаментальный компромисс между двумя классическими линейными детекторами в MIMO-системах.
Zero-Forcing полностью устраняет межпоточные интерференции, но ценой катастрофического усиления шума при плохой обусловленности канала. В нашем эксперименте при и SNR=8 дБ MSE детектора ZF составил 11.04 — созвездие превратилось в равномерное облако, структура QPSK полностью потеряна.
MMSE жертвует полным подавлением МСИ ради контроля над шумом. Регуляризатор предотвращает взрывной рост дисперсии в направлениях слабых собственных чисел канала. Результат: MSE = 0.50, что в 22 раза лучше ZF. На созвездиях это проявляется как чёткие компактные кластеры вокруг идеальных позиций QPSK.
Ключевой инсайт: преимущество MMSE над ZF максимально при низких SNR и плохой обусловленности канала — именно в тех условиях, где шум доминирует над сигналом. При высоких SNR или хорошей обусловленности оба детектора сходятся к псевдообращению , и их производительность практически идентична (что подтверждается совпадением BER-кривых в диапазоне 0-20 дБ для случайного Rayleigh-канала).
Практическая рекомендация: в реальных системах связи, где канал может временно становиться плохо обусловленным (замирания, корреляция антенн), MMSE предпочтительнее ZF как более робастный детектор. Вычислительная сложность обоих методов сопоставима (обе требуют обращения матрицы ), поэтому выбор в пользу MMSE не несёт дополнительных затрат.