精华内容
下载资源
问答
  • 正态分布参数估计
    千次阅读
    2020-10-25 18:46:36

    四、多元正态分布的参数估计

    1.多元正态分布的估计量

    对于多元正态分布 N p ( μ , Σ ) N_p(\mu,\Sigma) Np(μ,Σ),其参数只有两个——均值向量 μ \mu μ与自协方差矩阵 Σ \Sigma Σ,要对其进行估计,就要从总体中抽取简单随机样本。记抽取样本的容量为 n n n,每一个样本分别是 X ( α ) = ( x α 1 , ⋯   , x α p ) X_{(\alpha)}=(x_{\alpha1},\cdots,x_{\alpha p}) X(α)=(xα1,,xαp),将样本纵向排列,得到样本数据阵
    X = [ x 11 ⋯ x 1 p ⋮ ⋮ x n 1 ⋯ x n p ] . X=\begin{bmatrix} x_{11} & \cdots & x_{1p} \\ \vdots & & \vdots \\ x_{n1} & \cdots & x_{np} \end{bmatrix}. X=x11xn1x1pxnp.
    从样本数据阵出发,可以获得以下统计量:

    1. 样本均值 X ˉ \bar X Xˉ,这是对每个维度求均值,得到的一个 p p p维向量
      X ˉ = 1 n ∑ α = 1 n X ( α ) = ( x ˉ 1 , ⋯   , x ˉ p ) ′ = 1 n X ′ 1 n . \bar X=\frac 1n\sum_{\alpha=1}^n X_{(\alpha)}=(\bar x_1,\cdots ,\bar x_p)'=\frac 1nX'\boldsymbol 1_n. Xˉ=n1α=1nX(α)=(xˉ1,,xˉp)=n1X1n.
      这里 x ˉ i \bar x_i xˉi是对第 i i i个分量的平均,即
      x ˉ i = 1 n ∑ α = 1 n x α i . \bar x_i=\frac 1n\sum_{\alpha=1}^n x_{\alpha i}. xˉi=n1α=1nxαi.

    2. 样本离差阵 A A A,可以类比一维随机变量中的 ∑ i = 1 n ( x i − x ˉ ) 2 \sum_{i=1}^n (x_i-\bar x)^2 i=1n(xixˉ)2,即
      A = ∑ α = 1 n ( X ( α ) − X ˉ ) ( X ( α ) − X ˉ ) ′ A=\sum_{\alpha=1}^n(X_{(\alpha)}-\bar X)(X_{(\alpha)}-\bar X)' A=α=1n(X(α)Xˉ)(X(α)Xˉ)
      这样, A A A是一个 p × p p\times p p×p对角阵,它的第 ( i , j ) (i,j) (i,j)元,其实就是
      a i j = ∑ α = 1 n ( x α i − x ˉ i ) ( x α j − x ˉ j ) . a_{ij}=\sum_{\alpha=1}^n (x_{\alpha i}-\bar x_i)(x_{\alpha j}-\bar x_j). aij=α=1n(xαixˉi)(xαjxˉj).
      由此,还可以得到
      A = X ′ X − n X ˉ X ˉ ′ = X ′ [ I n − 1 n 1 n 1 n ′ ] X . A=X'X-n\bar X\bar X'=X'\left[I_n-\frac 1n\boldsymbol 1_n\boldsymbol 1_n' \right] X. A=XXnXˉXˉ=X[Inn11n1n]X.
      这个式子用来计算离差阵更为方便。

    3. 样本协方差阵 S S S,可以类比一维随机变量中的样本方差,即
      S = 1 n − 1 A , S=\frac 1{n-1}A, S=n11A,
      ( i , i ) (i,i) (i,i)元是变量 X i X_i Xi的样本方差,即
      s i i = 1 n − 1 ∑ α = 1 n ( x α i − x ˉ i ) 2 . s_{ii}=\frac 1{n-1}\sum_{\alpha=1}^n (x_{\alpha i}-\bar x_i)^2. sii=n11α=1n(xαixˉi)2.
      类似一维中样本方差的定义,也有
      S ∗ = 1 n ∑ α = 1 n ( x α i − x ˉ i ) 2 . S^*=\frac 1n\sum_{\alpha=1}^n(x_{\alpha i}-\bar x_i)^2. S=n1α=1n(xαixˉi)2.

    4. 样本相关阵 R R R,自然是由样本相关系数 r i j r_{ij} rij构成的 p × p p\times p p×p矩阵,即
      R = s i j s i i s j j = a i j a i i a j j . R=\frac{s_{ij}}{\sqrt{s_{ii}s_{jj}}}=\frac{a_{ij}}{\sqrt{a_{ii}a_{jj}}}. R=siisjj sij=aiiajj aij.

    有了这些统计量,我们就可以对总体的参数 μ , Σ \mu,\Sigma μ,Σ进行估计,使用的方法是最大似然估计。

    2.最大似然估计

    最大似然估计指的是,以使获得样本的出现几率最大的那组参数估计量,作为参数的点估计量。与一元情形类似,可以建立似然函数的概念。使用拉直运算,对 V e c ( X ′ ) {\rm Vec}(X') Vec(X)的密度函数建立似然函数,称为样本 X ( i ) X_{(i)} X(i)的似然函数(对数似然函数)。
    L ( μ , Σ ) = ∏ α = 1 n 1 ( 2 π ) p / 2 ∣ Σ ∣ 1 / 2 exp ⁡ [ − 1 2 ( x ( α ) − μ ) ′ Σ − 1 ( x ( α ) − μ ) ] = 1 ( 2 π ) n p / 2 ∣ Σ ∣ n / 2 exp ⁡ [ − 1 2 ∑ α = 1 n ( x ( α ) − μ ) ′ Σ − 1 ( x ( α ) − μ ) ] l ( μ , Σ ) = − n p 2 ln ⁡ ( 2 π ) + n 2 ln ⁡ ∣ Σ − 1 ∣ − 1 2 ∑ α = 1 n ( x ( α ) − μ ) ′ Σ − 1 ( x ( α ) − μ ) \begin{aligned} L(\mu,\Sigma)=&\prod_{\alpha=1}^n \frac{1}{(2\pi)^{p/2}|\Sigma|^{1/2}}\exp\left[-\frac12(x_{(\alpha)}-\mu)'\Sigma^{-1}(x_{(\alpha)}-\mu) \right] \\ =&\frac{1}{(2\pi)^{np/2}|\Sigma|^{n/2}}\exp\left[-\frac12\sum_{\alpha=1}^n(x_{(\alpha)}-\mu)'\Sigma^{-1}(x_{(\alpha)}-\mu) \right]\\ l(\mu,\Sigma)=&-\frac{np}2\ln(2\pi)+\frac n2\ln |\Sigma^{-1}|-\frac12\sum_{\alpha=1}^n(x_{(\alpha)}-\mu)'\Sigma^{-1}(x_{(\alpha)}-\mu) \end{aligned} L(μ,Σ)==l(μ,Σ)=α=1n(2π)p/2Σ1/21exp[21(x(α)μ)Σ1(x(α)μ)](2π)np/2Σn/21exp[21α=1n(x(α)μ)Σ1(x(α)μ)]2npln(2π)+2nlnΣ121α=1n(x(α)μ)Σ1(x(α)μ)

    要求其极大似然估计,需要对矩阵 Σ \Sigma Σ,向量 μ \mu μ求导(参见矩阵微商),得
    d l ( μ , Σ ) d μ = 1 2 ∑ α = 1 n ( Σ − 1 + ( Σ − 1 ) ′ ) ( x ( α ) − μ ) = Σ − 1 ( ∑ α = 1 n ( x ( α ) − μ ) ) = n Σ − 1 ( X ˉ − μ ) . d l ( μ , Σ ) d Σ − 1 = − n 2 Σ − 1 2 ∑ α = 1 n ( x ( α ) − μ ) ( x ( α ) − μ ) ′ = − 1 2 ( n Σ − A ) . \frac{{\rm d}l(\mu,\Sigma)}{{\rm d}\mu}=\frac12\sum_{\alpha=1}^n(\Sigma^{-1}+(\Sigma^{-1})')(x_{(\alpha)}-\mu)=\Sigma^{-1}(\sum_{\alpha=1}^n(x_{(\alpha)}-\mu))=n\Sigma^{-1}(\bar X-\mu).\\ \frac{{\rm d}l(\mu,\Sigma)}{{\rm d}\Sigma^{-1}}=-\frac n2\Sigma-\frac12\sum_{\alpha=1}^n(x_{(\alpha)}-\mu)(x_{(\alpha)}-\mu)'=-\frac12(n\Sigma-A). dμdl(μ,Σ)=21α=1n(Σ1+(Σ1))(x(α)μ)=Σ1(α=1n(x(α)μ))=nΣ1(Xˉμ).dΣ1dl(μ,Σ)=2nΣ21α=1n(x(α)μ)(x(α)μ)=21(nΣA).
    所以
    μ ^ = X ˉ , Σ ^ = A n . \hat \mu=\bar X,\quad \hat\Sigma = \frac An. μ^=Xˉ,Σ^=nA.

    用到的矩阵微商结论:对于对称阵 A A A与列向量 x x x,有
    d ln ⁡ ∣ A ∣ d A = A − 1 , d x ′ A x d A = x x ′ , d x ′ A x d x = ( A + A ′ ) x . \frac{{\rm d}\ln |A|}{{\rm d}A}=A^{-1},\\ \frac{{\rm d}x'Ax}{{\rm d}A}=xx',\\ \frac{{\rm d}x'Ax}{{\rm d}x}=(A+A')x. dAdlnA=A1,dAdxAx=xx,dxdxAx=(A+A)x.

    如果在已知 μ = μ 0 \mu=\mu_0 μ=μ0的情况下,依照以上过程,就可以得到
    Σ ^ = 1 n ∑ α = 1 n ( x ( α ) − μ 0 ) ( x ( α ) − μ 0 ) ′ . \hat \Sigma=\frac{1}{n}\sum_{\alpha=1}^n(x_{(\alpha)}-\mu_0)(x_{(\alpha)}-\mu_0)'. Σ^=n1α=1n(x(α)μ0)(x(α)μ0).
    所以,我们要找到 ( μ , Σ ) (\mu,\Sigma) (μ,Σ)的估计,就需要计算 ( X ˉ , A ) (\bar X,A) (Xˉ,A),接下来对它们进行性质讨论。

    3.最大似然估计的性质

    ( X ˉ , A ) (\bar X,A) (Xˉ,A)的分布具有类似一元统计中 X ˉ \bar X Xˉ S 2 S^2 S2的性质。

    定理:设 X ˉ \bar X Xˉ A A A分别是 p p p元正态总体 N p ( μ , Σ ) N_p(\mu,\Sigma) Np(μ,Σ)的样本均值向量和样本离差阵,则有

    1. X ˉ ∼ N p ( μ , 1 n Σ ) \bar X\sim N_p(\mu,\frac1n\Sigma) XˉNp(μ,n1Σ)
    2. A = d ∑ t = 1 n Z t Z t ′ A\stackrel {\rm d}=\sum\limits_{t=1}^n Z_tZ_t' A=dt=1nZtZt,其中 Z 1 , ⋯   , Z n − 1 Z_1,\cdots,Z_{n-1} Z1,,Zn1独立同 N p ( 0 , Σ ) N_p(0,\Sigma) Np(0,Σ)分布;
    3. X ˉ \bar X Xˉ A A A相互独立;
    4. P { A > 0 } = 1 ⇔ n > p {\rm P}\{A>0\}=1\Leftrightarrow n>p P{A>0}=1n>p

    前三个性质的证明方式也与一元情况类似,设 X X X为从多元正态总体中抽取的 n × p n\times p n×p样本数据阵, Γ \Gamma Γ n n n正交阵,形式如同
    Γ = [ r 11 ⋯ r 1 n ⋮ ⋮ r ( n − 1 ) 1 ⋯ r ( n − 1 ) n 1 / n ⋯ 1 / n ] = ( r i j ) n × n . \Gamma=\begin{bmatrix} r_{11} & \cdots & r_{1n} \\ \vdots & & \vdots \\ r_{(n-1)1} & \cdots & r_{(n-1)n} \\ 1/\sqrt n & \cdots & 1/\sqrt n \end{bmatrix}=(r_{ij})_{n\times n}. Γ=r11r(n1)11/n r1nr(n1)n1/n =(rij)n×n.

    Z = [ Z 1 ′ ⋮ Z n ′ ] = Γ [ X ( 1 ) ′ ⋮ X ( n ) ′ ] = Γ X . Z=\begin{bmatrix} Z_1' \\ \vdots \\ Z_n' \end{bmatrix} = \Gamma\begin{bmatrix} X_{(1)}' \\ \vdots \\ X_{(n)}' \end{bmatrix}=\Gamma X. Z=Z1Zn=ΓX(1)X(n)=ΓX.
    Z i ′ = ( r i 1 , ⋯ r i n ) X Z_i'=(r_{i1},\cdots r_{in})X Zi=(ri1,rin)X
    Z i = ( X ( 1 ) , ⋯   , X ( n ) ) [ r i 1 ⋮ r i n ] , i = 1 , ⋯   , n . Z_i=(X_{(1)},\cdots,X_{(n)})\begin{bmatrix} r_{i1} \\ \vdots \\ r_{in} \end{bmatrix},\quad i=1,\cdots,n. Zi=(X(1),,X(n))ri1rin,i=1,,n.
    因为 Z i Z_i Zi X ( 1 ) , ⋯   , X ( n ) X_{(1)},\cdots,X_{(n)} X(1),,X(n)的线性组合,所以 Z i Z_i Zi也是 p p p维正态向量,且
    E Z i = ∑ α = 1 n r i α E ( X ( α ) ) = { n ∑ α = 1 n r i α r n α μ = 0 , t ≠ n ; ∑ α = 1 n 1 n μ = n μ , t = n . C o v ( Z α , Z β ) = ∑ i = 1 n r α i r β i Σ = { O , α ≠ β ; Σ , α = β . {\rm E}Z_i=\sum_{\alpha=1}^n r_{i\alpha}{\rm E}(X_{(\alpha)})=\left\{ \begin{array}l \sqrt{n}\sum\limits_{\alpha=1}^n r_{i\alpha}r_{n\alpha}\mu=0,&t\ne n;\\ \sum\limits_{\alpha=1}^n \frac 1{\sqrt n}\mu=\sqrt n \mu,&t=n. \end{array} \right.\\ {\rm Cov}(Z_\alpha,Z_{\beta})=\sum_{i=1}^nr_{\alpha i}r_{\beta i}\Sigma=\left\{ \begin{array}l O,&\alpha\ne \beta;\\ \Sigma,&\alpha=\beta. \end{array} \right. EZi=α=1nriαE(X(α))=n α=1nriαrnαμ=0,α=1nn 1μ=n μ,t=n;t=n.Cov(Zα,Zβ)=i=1nrαirβiΣ={O,Σ,α=β;α=β.
    而显然 Z n = n X ˉ Z_n=\sqrt n\bar X Zn=n Xˉ,且 Z n ∼ N p ( n μ , Σ ) Z_n\sim N_p(\sqrt n\mu,\Sigma) ZnNp(n μ,Σ),所以 X ˉ ∼ N p ( μ , Σ / n ) \bar X\sim N_p(\mu,\Sigma/n) XˉNp(μ,Σ/n)。而
    ∑ α = 1 n Z α Z α ′ = ( Z 1 , ⋯   , Z n ) [ Z 1 ⋮ Z n ] = Z ′ Z = X ′ X , ∑ α = 1 n − 1 Z α Z α ′ = X ′ X − Z n Z n ′ = X ′ X − n X ˉ X ˉ ′ = A . \sum_{\alpha=1}^nZ_{\alpha}Z_{\alpha}'=(Z_1,\cdots,Z_n)\begin{bmatrix} Z_1\\ \vdots \\ Z_n \end{bmatrix}=Z'Z=X'X,\\ \sum_{\alpha=1}^{n-1}Z_{\alpha}Z_{\alpha}'=X'X-Z_nZ_n'=X'X-n\bar X\bar X'=A. α=1nZαZα=(Z1,,Zn)Z1Zn=ZZ=XX,α=1n1ZαZα=XXZnZn=XXnXˉXˉ=A.
    可以注意到, A A A Z 1 , ⋯   , Z n − 1 Z_1,\cdots,Z_{n-1} Z1,,Zn1的函数, X ˉ \bar X Xˉ Z n Z_n Zn的函数,又因为 Z 1 , ⋯   , Z n Z_1,\cdots,Z_n Z1,,Zn互相独立,所以 X ˉ \bar X Xˉ A A A相互独立。至于第四个性质,只需要记住,样本够多就能保证 A A A的非负定性即可。

    除此以外, X ˉ , A \bar X,A Xˉ,A作为 μ , Σ \mu,\Sigma μ,Σ的最大似然估计原型,还具有以下的性质:

    1. 无偏性: X ˉ \bar X Xˉ μ \mu μ的无偏估计, A / n A/n A/n不是 Σ \Sigma Σ的无偏估计,但 S = A / ( n − 1 ) S=A/(n-1) S=A/(n1) Σ \Sigma Σ的无偏估计。
    2. 有效性: X ˉ \bar X Xˉ S S S μ , Σ \mu,\Sigma μ,Σ的一致最小方差无偏估计,即 X ˉ , S \bar X,S Xˉ,S μ , Σ \mu,\Sigma μ,Σ的有效估计量。
    3. 相合性:当 n → ∞ n\to \infty n时, X ˉ , Σ ^ = A / n \bar X,\hat \Sigma=A/n Xˉ,Σ^=A/n μ , Σ \mu,\Sigma μ,Σ的强相合估计,即随着抽样数的增加,它们总会收敛于参数。
    4. 充分性: X ˉ , Σ ^ \bar X,\hat \Sigma Xˉ,Σ^ μ , Σ \mu,\Sigma μ,Σ的充分统计量。

    最大似然估计满足对参数函数依然适用的性质,即对于 μ , Σ \mu,\Sigma μ,Σ的最大似然估计 μ ^ , Σ ^ \hat \mu,\hat \Sigma μ^,Σ^,参数的函数 φ ( μ , Σ ) \varphi(\mu,\Sigma) φ(μ,Σ)的最大似然估计还是 φ ( μ ^ , Σ ^ ) \varphi(\hat \mu,\hat \Sigma) φ(μ^,Σ^)

    回顾总结

    1. 参数估计中,最重要的两个统计量是样本均值 X ˉ \bar X Xˉ与样本离差阵 A A A,它们与样本数据阵 X X X的关系分别是
      X ˉ = 1 n X ′ 1 n , A = X ′ X − n X ˉ X ˉ ′ = X ′ [ I n − 1 n 1 n 1 n ′ ] X . \bar X=\frac 1nX'\boldsymbol 1_n,\\ A=X'X-n\bar X\bar X'=X'\left[I_n-\frac1n\boldsymbol 1_n\boldsymbol 1_n' \right]X. Xˉ=n1X1n,A=XXnXˉXˉ=X[Inn11n1n]X.
      还有相关的统计量如 S = A / ( n − 1 ) , S ∗ = A / n S=A/(n-1),S^*=A/n S=A/(n1),S=A/n和样本相关阵 R , r i j = a i j / a i i a j j R,r_{ij}=a_{ij}/\sqrt{a_{ii}a_{jj}} R,rij=aij/aiiajj

    2. N p ( μ , Σ ) N_p(\mu,\Sigma) Np(μ,Σ)的参数 μ , Σ \mu,\Sigma μ,Σ的最大似然估计分别是 μ ^ = X ˉ , Σ ^ = A / n \hat \mu=\bar X,\hat \Sigma=A/n μ^=Xˉ,Σ^=A/n,一般可以由 X ˉ , A \bar X,A Xˉ,A估计出 μ , Σ \mu,\Sigma μ,Σ估计出。

    3. 关于 X ˉ , A \bar X,A Xˉ,A的性质,有 X ˉ , A \bar X,A Xˉ,A相互独立,且
      X ˉ ∼ N p ( μ , Σ / n ) , A = d ∑ α = 1 n − 1 Z α Z α ′ , Z α ∼ i . i . d . N p ( 0 , Σ ) . \bar X\sim N_p(\mu,\Sigma/n),\\ A\stackrel {\rm d}=\sum_{\alpha=1}^{n-1} Z_\alpha Z_{\alpha}',\quad Z_\alpha\stackrel {\rm i.i.d.}\sim N_p(0,\Sigma). XˉNp(μ,Σ/n),A=dα=1n1ZαZα,Zαi.i.d.Np(0,Σ).

    4. 在无偏性方面, X ˉ \bar X Xˉ μ \mu μ的无偏估计, S = A / ( n − 1 ) S=A/(n-1) S=A/(n1) Σ \Sigma Σ的无偏估计。

    5. 在有效性方面, X ˉ , S \bar X,S Xˉ,S μ , Σ \mu,\Sigma μ,Σ的最小方差无偏估计,即有效估计。

    6. 在相合性方面, X ˉ , A / n \bar X,A/n Xˉ,A/n μ , Σ \mu,\Sigma μ,Σ的强相合估计。

    7. 对于参数函数 φ ( μ , Σ ) \varphi(\mu,\Sigma) φ(μ,Σ),它的最大似然估计是 φ ( X ˉ , A / n ) \varphi(\bar X,A/n) φ(Xˉ,A/n)

    更多相关内容
  • 为了从观察图像中估计出合成孔径雷达(SAR)图像的雷达横截面积(RCS)模型的对数—正态分布的参数, 提出了强度图像和幅值图像的参数估计方法。推导出斑点的前两阶对数累积量和对数—正态分布的前两阶对数累积量; ...
  • 正态分布参数估计

    千次阅读 2020-12-30 21:26:12
    title: 正态分布参数估计 categories: 杂项 tags: 学习 正态分布参数估计 Untitled预处理clearrng(6331);mu = 1;sigma = 1;真实概率密度曲线:fplot(@(x) exp(-(x-mu).^2./(2*sigma))./(sqrt(2*pi)*sigma)...

    title: 正态分布与参数估计
    categories:

    • 杂项
      tags:
    • 学习

    正态分布与参数估计


    Untitled

    预处理

    clear
    rng(6331);
    mu = 1;
    sigma = 1;

    真实概率密度曲线:

    fplot(@(x) exp(-(x-mu).^2./(2*sigma))./(sqrt(2*pi)*sigma) )
    title("真实概率密度曲线")

    四个样本直方图

    n_1 = 10;
    n_2 = 50;
    n_3 = 500;
    n_4 = 1000;
    y_1 = sigma.*randn(n_1,1) + mu;
    y_2 = sigma.*randn(n_2,1) + mu;
    y_3 = sigma.*randn(n_3,1) + mu;
    y_4 = sigma.*randn(n_4,1) + mu;
    histogram(y_1,10)
    title("n=10")
    histogram(y_2,10)
    title("n=50")
    histogram(y_3,20)
    title("n=500")
    histogram(y_4,20)
    title("n=1000")

    四个样本参数估计值

    矩估计==最大似然估计
    [mu_, sigma_] =estimate(y_1 , n_1)
    mu_ = 0.6865
    sigma_ = 0.2972
    [mu_, sigma_] =estimate(y_2 , n_2)
    mu_ = 0.7853
    sigma_ = 1.0036
    [mu_, sigma_] =estimate(y_3 , n_3)
    mu_ = 0.9213
    sigma_ = 1.0244
    [mu_, sigma_] =estimate(y_4 , n_4)
    mu_ = 0.9908
    sigma_ = 1.0422

    参数估计误差随样本容量增大的变化曲线

    矩估计==最大似然估计
    样本容量取 10~1000
    n = 10:1:1000;
    dsigma = zeros(length(n),1);
    dmu = zeros(length(n),1);
    for i = 1:length(n)
    y = sigma.*randn(n(i),1) + mu;
    [mu_, sigma_] =estimate(y , n(i));
    dsigma(i) = abs(sigma - sigma_);
    dmu(i) = abs(mu - mu_);
    end
    plot(n,dsigma)
    hold on
    plot(n,dmu)
    legend('△σ','△μ')
    hold off

    矩估计、最大似然估计函数

    function [mu_, sigma_] = estimate(y , n)
    mu_ = sum(y)/n;
    sigma_ = sum((y-mu_).^2)/n;
    end

    展开全文
  • 第一章随机事件与概率 1第二章随机变量及概率分布 2第三章随机变量的数字特征 4第四章大数定律与中心极限定理 5第五章数理统计的基本概念 6第六章参数估计 8第
  • 在总体分布为N(μ,σ2)下,分别给出了未知参数μ,σ,σ2的几个无偏估计量,并加以证明。
  • 正态分布均值的贝叶斯估计公式的详细推导
  • 题目 1利用均匀分布U[-1,1],使用中心极限定理产生正态分布的随机数1000个,其中每次用来产生正态分布随机数的均匀分布样本容量个数 ≥ 50。 理论分析 由中心极限定理可知, 随机变量X1,X2,......Xn独立同分布,...

    本文包括以下内容:

    • 三道题的分析与代码实现
    • 实验值与理论值的比较
    • 完整代码附录

    题目 1  利用均匀分布U[-1,1],使用中心极限定理产生正态分布的随机数1000个,其中每次用来产生正态分布随机数的均匀分布样本容量个数 ≥ 50。

    理论分析


    由中心极限定理可知, 随机变量X1,X2, ......Xn独立同分布,且具有有限的数学期望和方差,则下面公式中的Yn将近似地服从标准正态分布N(0,1)。
     
    因此我们可以利用-1到1的均匀分布得到一系列Xi,随后分别产生1000个近似服从标准正态分布的单个随机数。

    MATLAB代码实现详解


    仅利用MATLAB中的rand()函数来产生均匀分布的随机数
    由于rand()函数仅产生(0,1)区间的均匀分布,因此还需做一些处理如下:
    function result = my_rand()
      below_zeros = -rand(1,50);
      above_zeros = rand(1,50);
      rand_sum = sum(below_zeros) + sum(above_zeros);
      total = 100;
      result = (rand_sum - total*0.5) / (10*total*total/12);
    end
    上方我们自定义的函数(非MATLAB内置的)是用于产生标准正态分布的单个随机数,所以循环执行该函数1000次并将其返回结果添加至数组之中就得到了我们需要的符合正态分布的随机数组rand_x。
    rand_x = [];
    for  time = 1:1000
        rand_x = [rand_x,my_rand()];
    end
     
    题目 2  在产生的1000个随机数中,均匀随机地抽取一个容量为100的样本。选择有效的点估计量给出基于该样本的总体均值和方差的点估计值,并将其与理论的推导的值比较。

    理论分析


    这一步涉及到了借助样本来估计总体的未知参数,也就是点估计。
    在这里采用了矩估计的方法,由于总体矩可以表现为未知参数的函数,我们用样本矩去替换总体矩,就构建了样本矩与未知参数的关系式。
    由以上公式可看出正态分布中的两个参数都可以由样本矩的式子表示出来,从而实现了对未知总体参数的估计。

    MATLAB代码实现详解


    第二题中需要随机抽取容量为100的样本。
    index = [];
    for time = 1:1000
        index = [index,time]; %产生1~1000的下标数组
    end
    rand_index = randperm(1000,100); %产生1~1000的100个不重复的随机下标
    sample = [];
    for time = 1:100
        sample = [sample,rand_x(rand_index(time))];
    end
    计算基于该样本估计的总体均值与方差:
    predict_ex = sum(sample)/100; %由样本估计出的均值
    sample1 = sample.^2; %将sample数组中的每一个值都平方
    square_sum = sum(sample1);
    predict_dx = square_sum/100 - predict_ex*predict_ex;
    %由样本估计的方差
     
    题目 3  基于抽取的容量的为100的样本给出总体均值和方差的置信度为0.95的置信区间,并考察总体参数是否落在置信区间内。

    理论分析


    这一步中要给出总体参数的置信区间。
    我们先在方差未知的情况下求均值的置信区间,有双侧置信区间公式如下图所示:
    该公式中n为样本容量100,由置信度=0.95,可知 α为0.05,查t分布表可知t0.025(99) = 1.9842
    然后在均值未知情况下求方差的置信区间, 有双侧置信区间公式如下图所示:
    该公式中n为样本容量100,由置信度=0.95,可知 α为0.05,同样的我们去查卡方分布表也能知道公式中有关量的具体值。

    MATLAB代码实现详解


    第三题中代入公式计算基于样本的总体均值和方差的置信区间。
    left_ex = predict_ex - 1.9842*sqrt(predict_dx)/10;
    %均值置信区间左侧
    right_ex = predict_ex + 1.9842*sqrt(predict_dx)/10;
    %均值置信区间右侧
    left_dx = (99*predict_dx)/128.422;
    %方差置信区间左侧
    right_dx = (99*predict_dx)/73.361;
    %方差置信区间右侧

    实验值与理论值比较


     
    第一次运行结果展示:
    可以看出由样本估计的均值为0.0103,由样本估计的方差为0.8851。
    均值的置信区间为(-0.1764,0.1970),方差的置信区间为(0.6823,1.1944)。
     
    第二次运行结果展示:
    可以看出由样本估计的均值为0.1013,由样本估计的方差为0.7502。
    均值的置信区间为(-0.0706,0.2731),方差的置信区间为(0.5783,1.0124)。
     
    总体均值与方差的理论值是0和1,这两次运行的样本估计值与理论值已经比较接近,且都落在了置信区间中。

    完整代码附录


     
    %产生标准正态分布,该部分涉及的my_rand()函数定义在最后
    rand_x = [];
    for time = 1:1000
        rand_x = [rand_x,my_rand()];
    end
    %样本抽取
    index = [];
    for time = 1:1000
        index = [index,time];
    end
    rand_index = randperm(1000,100);
    sample = [];
    for time = 1:100
        sample = [sample,rand_x(rand_index(time))];
    end
    %参数估计
    predict_ex = sum(sample)/100;
    sample1 = sample.^2;
    square_sum = sum(sample1);
    predict_dx = square_sum/100 - predict_ex*predict_ex;
    %求置信区间
    left_ex = predict_ex - 1.9842*sqrt(predict_dx)/10;
    right_ex = predict_ex + 1.9842*sqrt(predict_dx)/10;
    left_dx = (99*predict_dx)/128.422;
    right_dx = (99*predict_dx)/73.361;
    %展示结果数据
    disp(left_ex);
    disp(right_ex);
    disp(left_dx);
    disp(right_dx);
    disp(predict_ex);
    disp(predict_dx);
    %自定义函数用于产生单个正态分布随机数
    function result = my_rand()
      below_zeros = -rand(1,50);
      above_zeros = rand(1,50);
      rand_sum = sum(below_zeros) + sum(above_zeros);
      total = 100;
      result = (rand_sum - total*0) / (10*4/12);
    end
    展开全文
  • 文章讨论了正态分布中均值与方差变点问题,给出了变点位置的CUSUM型估计,证明了此估计是真实变点位置的强相合估计;与已有文献相比较,大大提高了收敛速度。
  • 对呈正态分布的样本的参数进行区间估计的方法,我们大致可以分成两类:去估计μ\muμ和去估计σ\sigmaσ,在估计μ\muμ 或σ\sigmaσ又分别有2种情况,则一共有以下四种情况: 已知σ\sigmaσ求μ\muμ 未知σ\...

    介绍

    在有参估计中,我们有两种常见的估计参数的方法:

    • 点估计
    • 区间估计

    本次估计方法,就是介绍一下用Python对正态分布样本进行区间估计的方法。

    数学方法概述

    对呈正态分布的样本的参数进行区间估计的方法,我们大致可以分成两类:去估计 μ \mu μ和去估计 σ \sigma σ,在估计 μ \mu μ σ \sigma σ又分别有2种情况,则一共有以下四种情况:

    1. 已知 σ \sigma σ μ \mu μ
    2. 未知 σ \sigma σ μ \mu μ
    3. 已知 μ \mu μ σ \sigma σ
    4. 未知 μ \mu μ σ \sigma σ

    面对4种不同的情况我们来分别简要讨论一下计算方法。

    1. 已知 σ \sigma σ μ \mu μ
      由于 σ \sigma σ已知,则我们设 D = X ‾ − μ σ / n ∼ N ( 0 , 1 ) D = \cfrac{\overline{X}-\mu}{\sigma/\sqrt{n}}\thicksim N(0,1) D=σ/n XμN(0,1),设某小概率事件为 α \alpha α,则置信度设置为 1 − α 1-\alpha 1α,为了尽可能缩小横坐标的取值范围,我们选择使用两边分别取 α 2 \cfrac {\alpha}{2} 2α的方法。
      这样 P { Φ − α 2 < D < Φ α 2 } = 1 − α P \{\Phi_{-\frac {\alpha}{2}} \lt D\lt \Phi_{\frac {\alpha}{2}} \} = 1-\alpha P{Φ2α<D<Φ2α}=1α
      化简求得,置信区间为:
      ( X ‾ − σ n Φ − α 2 , X ‾ + σ n Φ α 2 ) (\overline{X}-\cfrac {\sigma}{\sqrt{n}}\Phi_{-\frac {\alpha}{2}},\overline{X}+\cfrac {\sigma}{\sqrt{n}}\Phi_{\frac {\alpha}{2}}) (Xn σΦ2α,X+n σΦ2α)
      也就是说在这个区间里面取的值,由百分之 1 − α 1-\alpha 1α的可信度。
    2. 未知 σ \sigma σ μ \mu μ
      由于未知 σ \sigma σ,那么设 D = X ‾ − μ S / n ∼ t ( n − 1 ) D = \cfrac{\overline{X}-\mu}{S/\sqrt{n}}\thicksim t(n-1) D=S/n Xμt(n1)因为这里要进行无偏估计,t分布得自由度为 n − 1 n-1 n1,其中 S S S为样本方差。
      根据上面的方法同理求得,置信区间为: ( X ‾ − σ n t − α 2 , X ‾ + σ n t α 2 ) (\overline{X}-\cfrac {\sigma}{\sqrt{n}}t_{-\frac {\alpha}{2}},\overline{X}+\cfrac {\sigma}{\sqrt{n}}t_{\frac {\alpha}{2}}) (Xn σt2α,X+n σt2α)
    3. 已知 μ \mu μ σ \sigma σ
      已知 μ \mu μ,则设 D = ∑ i = 1 n ( x i − μ ) 2 σ 2 ∼ χ 2 ( n ) D = \cfrac{\sum_{i=1}^{n}(x_i-\mu)^2}{\sigma^2} \thicksim \chi^2(n) D=σ2i=1n(xiμ)2χ2(n)
      求解置信区间:
      ( ∑ i = 1 n ( x i − μ ) 2 χ 1 − α 2 2 , ∑ i = 1 n ( x i − μ ) 2 χ α 2 2 ) (\cfrac{\sum_{i=1}^{n}(x_i-\mu)^2}{\chi^2_{1-\frac{\alpha}{2}}},\cfrac{\sum_{i=1}^{n}(x_i-\mu)^2}{\chi^2_{\frac{\alpha}{2}}}) (χ12α2i=1n(xiμ)2,χ2α2i=1n(xiμ)2)
    4. 未知 μ \mu μ σ \sigma σ
      D = ( n − 1 ) S 2 σ 2 ∼ χ 2 ( n − 1 ) D = \cfrac{(n-1)S^2}{\sigma^2} \thicksim \chi^2(n-1) D=σ2(n1)S2χ2(n1),同样是为了保证无偏估计。
      求解置信区间为: ( ( n − 1 ) S 2 χ α 2 2 , ( n − 1 ) S 2 χ 1 − α 2 2 ) (\cfrac{(n-1)S^2}{\chi^2_{\frac{\alpha}{2}}},\cfrac{(n-1)S^2}{\chi^2_{1-\frac{\alpha}{2}}}) (χ2α2(n1)S2,χ12α2(n1)S2)

    估计期望的代码:

    import numpy as np
    import scipy.stats  as st
    from scipy.stats import t
    from scipy.stats import norm
    import matplotlib.pyplot as plt
    
    #在进行区间估计时,我们首先要知道样本总体的分布情况和一些参数
    #x表示数据取值的列表
    #sig表示方差
    #alpha一个足够小的概率
    #confidence表置信度
    def mean_ext(x,sig=None,confidence=0.95):
        alpha = 1 - confidence#求alpha
        if sig != None:#当方差已知时
            phi_alpha_2 = norm.isf(-alpha/2)#isf指的是残存函数的逆
            lower_limit = np.mean(x) - phi_alpha_2 * sig / np.sqrt(len(x)) #边界,最小值
            upper_limit = np.mean(x) + phi_alpha_2 * sig / np.sqrt(len(x))#边界,最大值
        else:#当方差未知时
            t_alpha_2 = t.isf(-alpha/2,df = len(x)-1)
            t_test = 1-t.isf(alpha/2,df = len(x)-1)
            std = np.std(x,ddof=1)
            lower_limit = np.mean(x) - t_alpha_2 * std / np.sqrt(len(x)-1)
            upper_limit = np.mean(x) + t_alpha_2 * std / np.sqrt(len(x)-1)
       
       
        return np.array(round(lower_limit, 6), round(upper_limit, 6))#返回一个四舍五入的值
    x = np.random.rand(500)生产数据
    xnew = np.linspace(0,10,10000) #画图数据
    interval_result = mean_ext(x,confidence=0.99) #区间估计结果,为一个区间
    u = np.mean(interval_result) #算个平均期望
    f = norm(loc = u,scale = np.std(x)).pdf(xnew) #pdf为概率密度函数
    f1 = norm(loc = u,scale = np.std(x)).pdf(x)
    plt.plot(x,f1,'r+',xnew,f,'g-')
    plt.legend(['point','pdf'],loc='best')
    plt.show()
    

    画图结果:在这里插入图片描述
    估计方差的代码:

    import numpy as np
    import scipy.stats  as st
    from scipy.stats import t
    from scipy.stats import norm
    import matplotlib.pyplot as plt
    
    #在进行区间估计时,我们首先要知道样本总体的分布情况和一些参数
    #x表示参数取值的列表
    #mean表期望
    #sig表示方差
    #alpha上分位
    def mean_ext(x,sig=None,confidence=0.95):
        alpha = 1 - confidence
        if sig != None:
            phi_alpha_2 = norm.isf(-alpha/2)
            lower_limit = np.mean(x) - phi_alpha_2 * sig / np.sqrt(len(x))
            upper_limit = np.mean(x) + phi_alpha_2 * sig / np.sqrt(len(x))
        else:
            t_alpha_2 = t.isf(alpha/2,df = len(x)-1)
            t_test = t.isf(1-alpha/2,df = len(x)-1)
            std = np.std(x,ddof=1)
            lower_limit = np.mean(x) - t_alpha_2 * std / np.sqrt(len(x)-1)
            upper_limit = np.mean(x) + t_alpha_2 * std / np.sqrt(len(x)-1)
       
       
        return np.array(round(lower_limit, 3), round(upper_limit, 3))
    
    
    def sigma_ext(x,mean = None,confidence = 0.99):
        alpha = 1-confidence
        top:float = 0#分子
        if mean != None:
            for xi in x:
                top += (xi-mean)**2
            chi2_alpha_2 = st.chi2.isf(1-alpha/2,df = len(x))
            chi2_alpha_2C = st.chi2.isf(alpha/2,df = len(x))
            chi2_alpha_test = chi2_alpha_2+chi2_alpha_2C
        else:
            d = np.mean(x)
            for xi in x:
                top +=(xi-d)**2
           
            chi2_alpha_2 = st.chi2.isf(1-alpha/2,df = len(x)-1)
            chi2_alpha_2C = st.chi2.isf(alpha/2,df = len(x)-1)
            chi2_alpha_test = chi2_alpha_2+chi2_alpha_2C
        lower_limit:float = top/chi2_alpha_2C
        upper_limit:float = top/chi2_alpha_2
    
        return np.array(round(lower_limit, 3), round(upper_limit, 3))
    
    #数据准备
    x = np.random.normal(0,1,300)
    xnew = np.linspace(-3,3,10000)
    interval_mean_result1 = mean_ext(x,confidence=0.99)
    interval_mean_result2 = mean_ext(x,1,confidence=0.99)
    interval_sigma_result1 = sigma_ext(x,confidence=0.99)
    interval_sigma_result2 = sigma_ext(x,0,confidence =0.99)
    u1 = np.mean(interval_mean_result1)
    u2 = np.mean(interval_mean_result2)
    sigma1 = np.mean(interval_sigma_result1)
    sigma2 = np.mean(interval_sigma_result2)
    f = norm(loc = 0,scale = 1).pdf(x)
    f1= norm(loc = u1,scale = np.std(x)).pdf(xnew)
    f2= norm(loc = u2,scale = 1).pdf(xnew)
    f3= norm(loc = np.mean(x),scale = np.sqrt(sigma1)).pdf(xnew)
    f4= norm(loc = 0,scale =np.sqrt(sigma2)).pdf(xnew)
    
    #画图模块
    plt.figure(figsize=(12,12))
    plt.subplot(221)
    plt.plot(x,f,'r+',xnew,f1,'g--')
    plt.legend(['point','pdf'],loc='best')
    plt.title("Mean Estimation(Unknown Sigma)")
    plt.subplot(222)
    plt.plot(x,f,'r+',xnew,f2,'g--')
    plt.legend(['point','pdf'],loc='best')
    plt.title("Mean Estimation(Known Sigma)")
    plt.subplot(223)
    plt.plot(x,f,'r+',xnew,f3,'b--')
    plt.legend(['point','pdf'],loc='best')
    plt.title("Sigma Estimation(Unknown Mean)")
    plt.subplot(224)
    plt.plot(x,f,'r+',xnew,f4,'b--')
    plt.legend(['point','pdf'],loc='best')
    plt.title("Sigma Estimation(Known Mean)")
    plt.show()
    
    

    图像:
    在这里插入图片描述

    展开全文
  • 函数 normfit()格式 [muhat,sigmahat,muci,sigmaci] = normfit(X) ; [muhat,sigmahat,muci,sigmaci] =normfit(X,alpha) 默认alpha=0.05即95%置信度 ...最大似然估计MLE(maximum likelihood estimatio...
  • 机器学习系列(二)多元正态分布

    千次阅读 2021-03-07 11:11:14
    一元正态分布回顾如果随机变量 服从均值为 方差为 的正态分布 (Univariate normal distribution), ,则其概率密度函数为:整个分布可以仅用均值及方差来刻画如果变量之间不相关,则它们相互独立经典统计检验通常...
  •  在实际应用中,多元正态分布中均值向量,和协差阵。通常是未知的,需由样本来估计,而参数估计方法很多,这里用最常见的最大似然估计法给出其估计量,并借助一元统计中学过的估计量性质指出这里给出的估计量也...
  • 正态分布时的贝叶斯估计

    千次阅读 2021-09-20 17:23:55
    正态分布时的贝叶斯估计 前导知识:【贝叶斯估计】 下面以最简单的一维正态分布模型为例来说明贝叶斯估计的应用。假设模型的均值μ\muμ是待估计参数,方差为σ2\sigma^2σ2为已知,分布密度写为: ρ(x∣μ)=12π...
  • 再求 ,使得似然函数最大]: 最大似然估计的性质: 由于充分利用了中的信息,极大似然估计 和 有着非常好的性质(估计的精确度很高),同时具有: 无偏性(需要修正),有效性,相合性,渐进正态性,且是 及 的...
  • 为了在统计过程中发现更多有趣的结果,我们将解决极大似然估计没有简单分析表达式的情况。举例来说,如果我们混合了各种分布
  • 混合正态分布用似然函数估计参数 图像为 数据为R中MASS包的geyser library(MASS) attach(geyser) #定义log-likelihood函数 ...
  • 2882.500000 50% 4250.000000 75% 5975.000000 max 24500.000000 Name: tradeMoney, dtype: float64 一个正态总体方差已知,均值的区间估计,使用的是正态分布 np.std 求得的均值是有偏的,这里我们需要的是无偏的...
  • 样本数据的类条件概率密度符合正态分布,对训练样本进行极大似然估计得到参数,再对测试样本进行分类。
  • 在做生物统计作业时用到了 Φ(-6.6),查表没有结果,所以想到用python解决。 以下代码主要是通过划分子区间求和的方式来计算。 import math ... return 0.5 #求标准正态分布的概率密度的积分 s=1/10000 x
  • 参数估计中,寿命数据是非常重要的. 传统的估计是基于完全精确的寿命数据. 然而,在实际中由于种种原因,有时收集的数据往往是不精确的. 这样,参数的模糊估计方法就十分必要. 本文将贝叶斯估计方法与模糊集理论相结合...
  • 许多作者已经基于非齐次集处理了正态分布参数估计:Hald A. 1949 [1],Arango-Castillo L.和Takahara G. 2018 [2]。 所有稳健的方法均基于这样的假设:受粗略误差影响的结果可以在检查点或截断点的左侧和/或右侧...
  • 多元正态分布参数估计PPT教案.pptx
  • 多元正态分布参数估计PPT学习教案.pptx
  • 多元正态分布参数估计与检验学习教案.pptx
  • 多元正态分布参数估计与检验PPT课件.pptx
  • EM算法在混合正态分布模型参数估计中的应用研究.docx
  • 多元正态分布参数估计与检验PPT学习教案.pptx
  • 正态分布简介你听说过钟形曲线吗?它往往是全球人们讨论最多的话题之一。很长一段时间以来,钟形曲线决定了对员工的专业评估,可以是一个受人喜爱或令人恐惧的话题,而这取决于与谁交谈!看看这张图片...
  • 正态分布的极大似然估计

    千次阅读 2022-02-28 13:12:47
    1. 正态分布的极大似然估计 笔记来源:Maximum Likelihood For the Normal Distribution, step-by-step!!! 1.1 正态分布参数对其形状的影响 1.1.1 μ值对正态分布的影响 1.1.2 σ值对正态分布的影响 1.2 极大...
  • 正态分布下的最大似然估计

    千次阅读 2021-09-19 15:11:18
    前导知识:【最大似然参数估计的求解】 本文仅以单变量正态分布情况下估计其均值和方差为例来说明最大似然估计的用法。 单变量正态分布的形式为: ρ(x∣θ)=12πσe−12(x−μσ)2(1) \rho(x|\theta)=\frac{1}{\...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 25,352
精华内容 10,140
关键字:

正态分布参数估计