随机产生高斯混合分布数据集

Halcom

matlab代码如下：

1. clc,clear,close all                            % 清屏、清工作区、关闭窗口
   
2. warning off                                 % 消除警告
   
3. feature jit off                               % 加速代码执行
   
4. % 产生X，并画图
   
5. randn('seed',0);                                 % 随机种子初始化
   
6. m1=[1, 2]'; m2=[2, 3]';m3=[3, 6]';                  % 已知量
   
7. m=[m1 m2 m3];
   
8. 
   
9. S(:,:,1)=0.1*eye(2);   % 初始化
   
10. S(:,:,2)=0.2*eye(2);   % 初始化
    
11. S(:,:,3)=0.3*eye(2);   % 初始化
    
12. P=[0.4 0.4 0.2];      % 初始化
    
13. N=1000;             % 初始化
    
14. sed=0;             % 初始化
    
15. [X,y]=mixt_model(m,S,P,N,sed);                  % 混合分布
    
16. plot_data(X,y,m,1)

复制代码

混合分布函数：

1. function [X,y]=mixt_model(m,S,P,N,sed)
   
2. % 产生一组数据，来自于混合高斯分布模型
   
3. %输入:
   
4. %   m:    均值向量
   
5. %   S:    协方差向量
   
6. %   P:    先验概率
   
7. %   N:    数据点数
   
8. %   sed:  rand()中的随机种子初始化
   
9. %输出:
   
10. %   X:  lxN 矩阵，列向量存储
    
11. %   y:  N为向量，对应于X下的正态分布指数
    
12. 
    
13. rand('seed',sed);   % 随机种子
    
14. [l,c]=size(m);      % 矩阵维数
    
15. P_acc=P(1);       % 初始化
    
16. for i=2:c
    
17.     t=P_acc(i-1)+P(i);
    
18.     P_acc=[P_acc t];
    
19. end
    
20. 
    
21. % 初始化，产生数据
    
22. X=[];
    
23. y=[];
    
24. for i=1:N
    
25.     t=rand;
    
26.     ind=sum(t>P_acc)+1;  % 正态分布指数
    
27.     X=[X; mvnrnd(m(:,ind)',S(:,:,ind),1)];
    
28.     y=[y ind];
    
29. end
    
30. X=X';

复制代码

绘图函数：

1. function plot_data(X,y,m,h)
   
2. % 画图，实现不同的数据类，用不同的颜色表示
   
3. % 输入:
   
4. %   X:  lxN matrix, 画图采用X的列向量进行一列一列绘制
   
5. %   y:  N-dimensional vector，数据类别
   
6. %   m:  lxc matrix, 均值向量.
   
7. %   h:  图像句柄
   
8. 
   
9. [l,N]=size(X); % N=no. of data vectors, l=dimensionality
   
10. [l,c]=size(m); % c=no. of classes
    
11. 
    
12. if(l~=2) || (c>7)
    
13.     fprintf('数据维数有问题，请检查\n')
    
14.     return
    
15. else
    
16.     pale=['r.'; 'g.'; 'b.'; 'y.'; 'm.'; 'c.';'co'];
    
17.     figure(h)
    
18.     % 画图
    
19.     hold on
    
20.     for i=1:N
    
21.         plot(X(1,i),X(2,i),pale(y(i),:))
    
22.         hold on
    
23.     end
    
24.    
    
25.     % 绘制数据类重心
    
26.     for j=1:c
    
27.         plot(m(1,j),m(2,j),'k+')
    
28.         hold on
    
29.     end
    
30. end

复制代码