MATLAB矩阵：rand/randn/randi 三兄弟：均匀分布、正态分布、随机整数的区别与用法

一、随机性的工程起点：为什么PRNG不是“真随机”

要理解rand、randn、randi的区别，必须先回到一个被多数教程跳过的前提：计算机无法凭空产生随机性。MATLAB提供的这三个函数，本质上都是伪随机数生成器（Pseudo-Random Number Generator, PRNG）的输出接口。所谓“伪”，是指它们由一个确定性算法在有限状态空间上迭代产生，给定相同种子，输出序列完全可复现。

MATLAB自R2007b起，默认PRNG切换为梅森旋转算法（Mersenne Twister, MT19937）。该算法由Matsumoto与Nishimura于1998年提出，周期为2¹⁹⁹³⁷−1，状态向量长度为624个32位字。这个周期长度意味着：即使每秒消耗10⁹个随机数，也需要约10⁶⁰⁰⁰年才会遍历完整个序列。本文评述：周期长度常被误读为“随机性好”，但周期只是必要条件而非充分条件——一个周期极长的序列仍可能通过不了谱检验或k维均匀性检验。

MT19937的初始化过程值得单独说明。当调用rng(seed)时，MATLAB并非简单地把种子直接填入状态数组，而是通过一个基于线性同余的初始化循环，将单个整数扩展为624个字的状态向量。这意味着：相邻种子（如1和2）产生的序列在统计上并不“相邻”，但也不保证完全独立。文献[3]指出，MT19937在初始化阶段若种子熵不足，前若干个输出可能存在弱相关性，这在密码学场景中不可接受，但在一般科学计算中影响可忽略。

那么，rand、randn、randi三者的关系是什么？简言之：它们共享同一个底层PRNG引擎，只是对原始32位整数输出做了不同的后处理变换。rand将整数映射到[0,1)区间的双精度浮点数；randn在均匀分布基础上施加变换以逼近标准正态分布；randi则将输出映射到指定整数区间。理解这一点，就能明白为什么切换rng种子会同时影响三者——因为它们消费的是同一条随机流。

拓展阅读：MathWorks官方文档《Random Number Generation》系统介绍了MATLAB的PRNG架构，地址：https://www.mathworks.com/help/matlab/random-number-generation.html 。此外，Melissa O'Neill的PCG随机数生成器主页（https://www.pcg-random.org/）提供了对MT19937局限性的深入讨论，适合进阶读者。

二、rand：均匀分布的算法内核与参数语义

2.1 数学定义与输出区间

rand的数学语义是产生服从连续均匀分布U(0,1)的样本。严格来说，MATLAB文档明确说明rand返回的是开区间(0,1)内的值，而非闭区间[0,1]。这一细节在涉及对数变换或除法运算时至关重要——若rand可能返回0，则log(rand())会返回-Inf，进而污染整个计算链。

从实现角度看，rand将MT19937产生的32位无符号整数除以2³²，再经过舍入处理映射到双精度浮点数的有效位数。由于双精度浮点数在[0.5,1)区间有2⁵²个可表示值，而在[0,0.5)区间同样有2⁵²个，因此rand的输出在浮点意义上并非严格均匀——极小值附近的密度略低于理论值。本文评述：这种偏差在绝大多数工程应用中可忽略，但在极端尾部概率估计（如p<10^-15）场景中需要警惕。

2.2 参数语义与常用调用形式

rand的调用形式看似简单，但参数语义存在易混淆之处：

• rand(n)：返回n×n矩阵，元素服从U(0,1)。
• rand(m,n)：返回m×n矩阵。
• rand(sz1,sz2,...,szN)：返回多维数组。
• rand('like',p)：返回与p同类型、同分布的数组（R2016b引入）。

需要特别注意的是，rand('seed',x)这类旧式语法在R2007b之后已被弃用，但仍能在部分遗留代码中见到。本文评述：旧语法不仅影响可读性，还会与rng机制产生冲突，建议在代码审查中将其列为“必须重构项”。

2.3 从U(0,1)到任意区间：线性变换的正确姿势

若需要产生[a,b]区间上的均匀分布，标准做法是a + (b-a)*rand()。但这里存在一个常被忽视的边界问题：由于rand返回(0,1)开区间，变换后无法精确取到端点a和b。在需要闭区间的场景（如随机索引），应改用randi。

% 产生[-5, 5]区间均匀分布，100万个样本
rng(42);  % 固定种子保证可复现
n = 1e6;
x = -5 + 10*rand(n,1);
fprintf('样本均值: %.4f (理论值 0)\n', mean(x));
fprintf('样本方差: %.4f (理论值 %.4f)\n', var(x), 100/12);
fprintf('最小值: %.6f, 最大值: %.6f\n', min(x), max(x));

运行上述代码（MATLAB R2024b，模拟数据），可观察到样本均值约在0附近波动，方差接近8.3333。这验证了线性变换保持均匀分布性质。但笔者提醒：当变换涉及非线性函数时，均匀性不再保持，这是下一节randn要解决的核心问题。

三、randn：正态分布背后的变换方法与尾部风险

3.1 正态分布为何不能“直接生成”

与均匀分布不同，正态分布没有有界支撑集，其概率密度函数在实数轴上非零。这意味着无法通过简单的线性映射从U(0,1)得到N(0,1)。MATLAB采用的是一种经典变换方法：Ziggurat算法（R2007b之后）或更早的Box-Muller变换。

Box-Muller变换的数学形式为：若U₁、U₂独立服从U(0,1)，则Z₁=√(−2lnU₁)·cos(2πU₂)、Z₂=√(−2lnU₁)·sin(2πU₂)独立服从N(0,1)。该变换涉及对数、三角函数和平方根，计算成本较高。Ziggurat算法则通过分层拒绝采样，将大部分样本的生成简化为查表和比较，速度提升约3-5倍。文献[7]的基准测试表明，在Intel Xeon平台上，Ziggurat生成10⁷个正态样本约耗时0.12秒，而Box-Muller约0.45秒。

本文评述：变换方法的选择不仅影响速度，还影响尾部精度。Box-Muller在极端尾部（|z|>6）可能因浮点舍入产生偏差，而Ziggurat通过精心设计的层边界，在尾部区域仍能保持较好的精度。对于金融风险建模等关注尾部事件的场景，这一差异不可忽视。

3.2 参数语义与常见误用

randn的调用形式与rand类似，但语义完全不同：

调用形式	输出分布	典型用途
randn(n)	N(0,1)，n×n矩阵	标准正态样本
randn(m,n)	N(0,1)，m×n矩阵	权重初始化
mu + sigma*randn(...)	N(mu, sigma²)	任意正态分布
randn('like',p)	与p同类型	GPU/分布式数组

最常见的误用是混淆方差与标准差。randn产生的是标准差为1、方差为1的标准正态分布。若需要方差为σ²的分布，应使用sigma*randn()而非sigma^2*randn()。笔者在代码审查中多次见到这一错误，其后果在深度学习初始化中尤为严重——方差被错误放缩会导致梯度爆炸或消失。

3.3 正态性检验：你的randn真的正态吗？

虽然MATLAB的randn经过严格测试，但在特定应用场景中，验证其正态性仍是良好实践。常用的检验方法包括：

• Shapiro-Wilk检验：适用于小样本（n<5000），对尾部偏差敏感。
• Kolmogorov-Smirnov检验：适用于大样本，检验经验分布与理论分布的最大偏差。
• Anderson-Darling检验：对尾部赋予更高权重，适合金融应用。
• 偏度与峰度：快速诊断，正态分布偏度为0、峰度为3。

% 正态性检验示例（模拟数据）
rng(2024);
x = randn(10000,1);

% 偏度与峰度
skewness_val = skewness(x);
kurtosis_val = kurtosis(x);
fprintf('偏度: %.4f (理论值 0)\n', skewness_val);
fprintf('峰度: %.4f (理论值 3)\n', kurtosis_val);

% Kolmogorov-Smirnov检验
[h,p] = kstest((x-mean(x))/std(x), 'CDF', makedist('Normal'));
fprintf('KS检验 p值: %.4f (p>0.05则不能拒绝正态假设)\n', p);

运行结果（模拟数据）通常显示偏度接近0、峰度接近3，KS检验p值大于0.05。但笔者提醒：通过正态性检验不等于“完美正态”，检验的效力取决于样本量和显著性水平，大样本下微小偏差也会导致拒绝。

四、randi：随机整数的区间语义与常见误用

4.1 离散均匀分布的定义

randi产生的是离散均匀分布的随机整数，即区间[1, imax]或[imin, imax]内每个整数以等概率出现。其数学定义为P(X=k)=1/(imax−imin+1)，k∈{imin,...,imax}。

与rand的连续均匀分布不同，randi的输出是离散的，且包含两个端点。这一特性使其成为随机索引、随机抽样、随机分组等任务的首选工具。

4.2 参数语义详解

调用形式	输出范围	说明
randi(imax)	[1, imax]	标量，默认从1开始
randi(imax,n)	[1, imax]	n×n矩阵
randi([imin,imax],m,n)	[imin, imax]	指定区间
randi(imax,sz,'like',p)	[1, imax]	指定类型

一个常见误用是randi([0, n])与randi(n+1)-1的等价性。两者确实等价，但前者更直观、更不易出错。本文评述：在代码可读性优先的原则下，应优先使用显式区间形式。

4.3 随机抽样与随机分组实战

randi在工程中最典型的应用是随机抽样和随机分组。以下代码演示如何从1000个样本中随机抽取100个不重复样本：

% 无放回随机抽样（模拟数据）
rng(7);
N = 1000;  % 总体大小
n = 100;   % 样本量
idx = randperm(N, n);  % 推荐方法：randperm
% 等价但效率较低的方法：循环拒绝采样
% 不推荐：randi可能产生重复索引

fprintf('抽样索引范围: [%d, %d]\n', min(idx), max(idx));
fprintf('是否存在重复: %d\n', length(unique(idx)) ~= n);

笔者强调：无放回抽样应使用randperm而非randi，因为randi允许重复，需要额外的去重逻辑，既低效又易错。randperm内部同样基于PRNG，但保证了输出的唯一性。

五、三兄弟横向对比：一张表看清选型边界

经过前四节的拆解，我们可以从多个维度对三兄弟进行系统对比。下表是本文的核心总结之一：

维度	rand	randn	randi
分布类型	连续均匀U(0,1)	标准正态N(0,1)	离散均匀
输出范围	(0,1)开区间	(−∞,+∞)	[imin,imax]闭区间
输出类型	double	double	double（可指定整数类）
底层算法	MT19937+缩放	MT19937+Ziggurat	MT19937+取模/缩放
典型用途	概率抽样、蒙特卡洛	噪声建模、权重初始化	随机索引、分组
常见误用	当作闭区间使用	方差/标准差混淆	无放回抽样

本文评述：三者的选型不应基于“习惯”或“顺手”，而应基于目标分布的数学定义。若目标分布是连续均匀，选rand；若目标分布是正态，选randn；若目标分布是离散均匀，选randi。任何“用rand模拟randn”或“用randi模拟rand”的做法，都是在引入不必要的统计偏差。

六、统计检验：如何验证你的随机数“合格”

6.1 均匀性检验：卡方与KS

对于rand产生的均匀分布，最直接的检验是卡方拟合优度检验：将[0,1]划分为k个等宽区间，统计各区间频数，与理论频数比较。另一个常用方法是Kolmogorov-Smirnov检验，直接比较经验累积分布函数与理论CDF的最大偏差。

% 均匀性检验（模拟数据）
rng(123);
n = 1e5;
x = rand(n,1);

% 卡方检验：10个等宽区间
edges = 0:0.1:1;
counts = histcounts(x, edges);
expected = n/10 * ones(1,10);
chi2_stat = sum((counts - expected).^2 ./ expected);
p_chi2 = 1 - chi2cdf(chi2_stat, 9);
fprintf('卡方统计量: %.2f, p值: %.4f\n', chi2_stat, p_chi2);

% KS检验
[h,p_ks] = kstest(x, 'CDF', makedist('Uniform',0,1));
fprintf('KS检验 p值: %.4f\n', p_ks);

6.2 独立性检验：自相关与游程

均匀性只是随机性的一个维度，另一个关键维度是独立性。若相邻样本存在相关性，即使边际分布均匀，序列仍不是“好”的随机数。常用检验包括：

• 自相关函数（ACF）：计算滞后k的自相关系数，理想情况下应接近0。
• 游程检验（Runs Test）：检验序列中高于/低于中位数的连续段长度是否随机。
• 谱检验（Spectral Test）：检测高维格点结构，MT19937在维度较高时可能出现格点规律。

本文评述：对于一般工程应用，通过均匀性和低阶自相关检验即可满足需求。但对于蒙特卡洛积分等高精度场景，建议使用MATLAB的runstest和autocorr进行快速诊断。

七、可复现性：rng种子机制的深层逻辑

7.1 rng函数的三要素

MATLAB的rng函数是控制随机数生成的核心接口，其完整语义包含三个要素：种子（seed）、生成器类型（generator）、状态（state）。

% rng的典型用法
rng(42);                    % 设置种子为42，使用默认生成器
rng('default');             % 恢复默认设置（种子0，MT19937）
rng('shuffle');             % 基于当前时间设置种子
s = rng;                    % 保存当前状态
rng(s);                     % 恢复保存的状态

% 指定生成器类型
rng(42, 'twister');         % MT19937
rng(42, 'simdTwister');     % SIMD加速版MT19937
rng(42, 'combRecursive');   % 组合递归生成器
rng(42, 'philox');          % Philox 4x32（GPU推荐）
rng(42, 'threefry');        % Threefry 4x64

7.2 可复现性的边界条件

许多用户认为“设置rng种子就能保证结果完全可复现”，但这一认知存在多个边界条件：

• MATLAB版本差异：不同版本的MT19937实现可能有细微差异，跨版本复现需谨慎。
• 并行计算：parfor中的随机流分配与串行不同，需使用parallel.gpu.rng或parfor的随机流控制。
• GPU计算：gpuArray上的rand/randn使用不同的生成器，需单独设置。
• 函数内部调用：某些工具箱函数内部会调用rand，可能“消耗”随机流。

本文评述：可复现性不是“设置种子”这一动作，而是一个系统工程问题。在关键实验中，建议同时保存代码版本、MATLAB版本、rng状态和硬件信息，形成完整的“随机性审计链”。

八、工程场景实战：从蒙特卡洛到深度学习初始化

8.1 蒙特卡洛积分：rand的主场

蒙特卡洛积分是rand最经典的应用场景。以估计π为例，在单位正方形内随机投点，统计落在四分之一圆内的比例：

% 蒙特卡洛估计π（模拟数据）
rng(2024);
n = 1e7;
x = rand(n,1);
y = rand(n,1);
inside = sum(x.^2 + y.^2 <= 1);
pi_est = 4 * inside / n;
fprintf('π估计值: %.6f (真实值: 3.141593)\n', pi_est);
fprintf('相对误差: %.4f%%\n', abs(pi_est-pi)/pi*100);

运行结果（模拟数据）通常显示π估计值在3.1416附近，相对误差约0.01%。蒙特卡洛方法的误差收敛速度为O(1/√n)，这意味着要将误差降低一个数量级，需要增加100倍样本量。本文评述：蒙特卡洛的精度瓶颈不在rand的质量，而在样本量的平方根律，这是概率论的基本限制。

8.2 深度学习权重初始化：randn的精细调参

在深度学习领域，权重初始化对训练成败至关重要。Xavier初始化（Glorot & Bengio, 2010）和He初始化（He et al., 2015）都基于randn：

% Xavier初始化与He初始化（模拟数据）
rng(42);
fan_in = 512;
fan_out = 256;

% Xavier初始化：适用于tanh/sigmoid
W_xavier = randn(fan_in, fan_out) * sqrt(2/(fan_in+fan_out));
fprintf('Xavier: 标准差=%.4f, 理论值=%.4f\n', ...
    std(W_xavier(:)), sqrt(2/(fan_in+fan_out)));

% He初始化：适用于ReLU
W_he = randn(fan_in, fan_out) * sqrt(2/fan_in);
fprintf('He: 标准差=%.4f, 理论值=%.4f\n', ...
    std(W_he(:)), sqrt(2/fan_in));

本文评述：初始化方法的选择与激活函数强相关。ReLU的半区死亡特性要求更大的初始方差（He初始化），而tanh的饱和特性则要求更小的方差（Xavier初始化）。误用初始化方法是深度学习训练不收敛的常见原因之一。

8.3 通信系统仿真：randn建模AWGN信道

在通信系统仿真中，加性高斯白噪声（AWGN）信道是基本模型。randn用于生成噪声样本：

% AWGN信道误码率仿真（模拟数据）
rng(2024);
N = 1e6;           % 比特数
SNR_dB = 0:2:10;   % 信噪比范围
BER = zeros(size(SNR_dB));

for k = 1:length(SNR_dB)
    snr = 10^(SNR_dB(k)/10);
    bits = randi([0,1], N, 1);           % 随机比特
    symbols = 2*bits - 1;                % BPSK调制
    noise = randn(N,1) / sqrt(2*snr);    % AWGN噪声
    received = symbols + noise;
    decisions = received > 0;
    BER(k) = mean(decisions ~= bits);
end

disp(table(SNR_dB', BER', 'VariableNames', {'SNR_dB','BER'}));

运行结果（模拟数据）应与理论BPSK误码率Q(√(2·SNR))吻合。本文评述：randn的尾部精度直接影响高信噪比下的误码率估计，因为低误码率事件依赖于噪声分布的极端尾部。

九、并行与GPU：随机流的独立性与可重复性挑战

9.1 parfor中的随机数陷阱

在并行计算中，多个worker同时调用rand会产生随机流冲突——不同worker可能产生相同的随机序列，破坏统计独立性。MATLAB的解决方案是为每个worker分配独立的子流：

% parfor中的随机流控制（模拟数据）
rng(42);  % 主线程设置种子
n = 100;
results = zeros(n,1);

parfor i = 1:n
    % 每个worker基于i生成独立子流
    rng(42 + i, 'twister');
    results(i) = randn();
end

fprintf('并行结果均值: %.4f\n', mean(results));

本文评述：parfor中的随机数必须显式控制，否则结果不可复现且可能统计相关。MATLAB R2019b之后提供了更优雅的parallel.pool.Constant和RandStream机制。

9.2 GPU随机数生成

GPU上的随机数生成面临独特挑战：数千个线程需要同时获取随机数，且不能有冲突。MATLAB的Parallel Computing Toolbox为GPU提供了专门的生成器（如Philox和Threefry），这些生成器支持基于计数器的并行流，每个线程根据其全局ID生成独立序列。

% GPU随机数生成（需Parallel Computing Toolbox）
if canUseGPU
    gpu_rng = gpurng;              % 获取GPU随机流
    gpurng(42);                    % 设置GPU种子
    x_gpu = rand(1000, 'gpuArray'); % GPU上生成随机数
    x_cpu = gather(x_gpu);         % 传回CPU
    fprintf('GPU随机数均值: %.4f\n', mean(x_cpu));
end

文献[15]的研究表明，GPU随机数生成器的吞吐量可达CPU的50-100倍，但统计质量需要单独验证——某些早期GPU生成器在低位比特上存在规律性。

十、前沿演进与学术预判：从MT到PCG、从CPU到量子

10.1 MT19937的局限与替代方案

尽管MT19937是MATLAB的默认生成器，但学术界对其局限性的讨论从未停止。主要问题包括：

• 状态空间过大：624个32位字的状态向量在内存受限场景中不友好。
• 初始化成本高：每次rng调用需初始化624个字。
• 线性结构：MT19937是线性生成器，存在可预测性，不适合密码学。
• 高维格点：在维度超过623时，输出可能出现格点规律。

PCG（Permuted Congruential Generator）系列由Melissa O'Neill于2014年提出，结合了线性同余生成器的速度和统计质量更好的输出函数。PCG的周期可达2¹²⁸，状态仅需128位，且通过了更严格的统计测试。本文评述：PCG在状态大小、速度和统计质量之间取得了更好的平衡，是MT19937的有力替代者。MATLAB R2023a之后已支持通过第三方工具箱使用PCG。

10.2 量子随机数生成：真正的随机性？

量子随机数生成器（QRNG）利用量子力学的不确定性原理产生真随机数。与PRNG不同，QRNG的输出在理论上不可预测。近年来，基于光子、原子和超导电路的QRNG已实现商业化，如ID Quantique的Quantis系列。

然而，QRNG在科学计算中的普及面临成本、速度和接口兼容性等挑战。本文评述：对于绝大多数工程应用，高质量PRNG已足够；QRNG的价值主要体现在密码学、安全通信和量子模拟等对“真随机性”有严格要求的领域。

10.3 可微分随机性：深度学习的新需求

随着变分自编码器（VAE）和扩散模型的发展，可微分随机数生成成为新的研究热点。传统randn生成的样本不可微分，无法通过反向传播优化。重参数化技巧（Reparameterization Trick）通过z = μ + σ·ε（ε~N(0,1)）将随机性从计算图中分离，使梯度可以流过μ和σ。

文献[22]提出的Gumbel-Softmax技巧则将离散随机采样松弛为连续可微形式，广泛应用于离散隐变量模型。本文评述：可微分随机性正在重新定义“随机数”在优化中的角色——从不可控的噪声源转变为可训练的参数。

十一、总结：一条主线串起三兄弟的选型决策树

回顾全文，我们可以将rand、randn、randi的选型逻辑凝练为一棵决策树：

本文的核心论点可以概括为：随机函数的选择不是语法问题，而是统计建模问题。rand、randn、randi分别对应三种不同的概率分布假设，误用它们等价于在模型中引入了错误的先验。笔者建议读者在每次使用随机函数时，都问自己三个问题：目标分布是什么？参数是否正确？结果是否可复现？

最后，随机数生成是一个快速演进的领域。从MT19937到PCG，从CPU到GPU，从伪随机到量子随机，工具在变，但对统计假设的敬畏和对可复现性的追求始终是工程实践的基石。

十二、参考文献与拓展资源