（20）MATLAB使用卡方（chi-square）分布生成Nakagami-m分布

---

前言

Nakagami分布随机变量的生成，可以使用gamma分布实现，也可以使用卡方（chi-square）分布随机变量生成。使用伽马（gamma）分布生成Nakagami-m分布的方法见：
（17）MATLAB使用伽马（gamma）分布生成Nakagami-m分布的方法1
（18）MATLAB使用伽马（gamma）分布生成Nakagami-m分布的方法2

本文介绍使用卡方（chi-square）分布随机变量生成Nakagami分布随机变量的方法。

---

一、使用卡方分布函数生成Nakagami分布随机变量

使用卡方分布生成Nakagami分布随机变量的方法由下式给出：

为了得到给定形状参数（m）和扩展参数（ω）的Nakagami-m随机变量，可将MATLAB仿真分为如下步骤：
（1）首先生成卡方分布随机变量CHI2，其自由度k设置为2m；
（2）将生成的随机变量CHI2乘以ω/（2m）并开方。
下面给出使用该方法生成Nakagami分布随机变量的MATLAB代码。

二、MATLAB仿真

以下代码使用卡方分布生成不同的形状参数（m）和尺度参数（ω）的Nakagami随机变量。然后给出仿真Nakagami随机变量的PDF画图结果。

1.仿真代码

代码如下：

clc
clear all
close all

m = 1;                         % Nakagami分布的形状参数（shape parameters）
Omega = 0.2;                   % Nakagami分布的尺度参数（spread parameters）
N = 1e7;                       % Number of Samples

% 生成卡方分布随机数
% Generate Chi-squared distributed random numbers
% 方法1：使用标准正态分布随机变量的平方和生成卡方分布随机数
CHI2 = 0;                       % 卡方分布随机变量（Chi-squared distributed random numbers）
k = 2*m;                        % 卡方分布的自由度（degrees of freedom）
for j=1:k
    % Sum of square of k independent standard normal Random Variables 
    % k个独立标准正态随机变量的平方和
    CHI2 = CHI2 + randn(1 ,N).^2;
end

% 方法2：使用库函数chi2rnd直接生成卡方分布随机数
% CHI2 =  chi2rnd(2*m,1,N);

% 使用卡方分布随机变量生成Nakagami分布的随机变量
Y = sqrt(Omega/(2*m)) * sqrt(CHI2);

% 画图方法1。概率密度函数（使用了梯形积分）
figure();
[Q,X] = hist(Y,1000);           % histogram，Q为划分到每个桶中的元素个数，X为每个区间的中心
plot(X,Q/trapz(X,Q),'r');       % plot estimated PDF，trapz(X,Q)使用梯形法计算Q相对于X的积分
grid on;
title('Nakagami-m PDF ' );
xlabel('Parameter - y' );       % y应该是随机变量
ylabel('f_Y(y)' );              % 随机变量的概率密度

% 画图方法2。绘制生成的随机数的直方图，以概率密度形式（MATLAB建议的方式）
nbins = 50;                     % bin数量
figure();
histogram(Y,nbins,'Normalization','pdf','DisplayStyle','bar');
hold on

%% nakagami分布的概率密度的理论值
% Nakagami分布的参数
% m = 1;
% Omega = 2;

% nakagami分布的概率密度的理论表达式
x = 0.01:0.01:3;
f = (2*m.^m./(gamma(m).*Omega^m)) .* x.^(2*m-1) .* exp(-m*x.^2./Omega);
% figure(1)
plot(x,f,'LineWidth',1.5)
grid on
title('Nakagami分布的概率密度')
legend('概率密度函数的估计值','概率密度函数的理论值')

2.运行结果

代码运行结果画图如下：

---