如何剔除MATLAB高斯Copula轮廓线外的散点？

现有MATLAB代码基于高斯Copula函数生成包含轮廓线内外散点的图形，其中givenData是5000×2的矩阵。需要修改代码，仅绘制最外层红色轮廓线内的散点，剔除轮廓外的点。

当前绘图代码片段 #

plot(givenData(:,1),givenData(:,2),'b.','MarkerSize',3);

双变量轮廓图函数原代码 #

function [xgrid,ygrid,Z] = biVariateContourPlotsGMMCopula(givenData,gmmObject,~,numMeshPoints,x_dim,y_dim)

d = 2;
if nargin < 5
    x_dim = 1;
    y_dim = 2;
end

if x_dim == y_dim
    hist(givenData(:,x_dim),10);
    return;
end

numMeshPoints = min(numMeshPoints,256);

givenData = givenData(:,[x_dim y_dim]);
alpha = gmmObject.alpha;
mu = gmmObject.mu(:,[x_dim y_dim]);
sigma = gmmObject.sigma([x_dim y_dim],[x_dim y_dim],:) + 0.005*repmat(eye(d),[1 1 numel(alpha)]);


gmmObject = gmdistribution(mu,sigma,alpha);

bin_num = 256;
for j = 1:2
   l_limit = min(gmmObject.mu(:,j))-3*(max(gmmObject.Sigma(j,j,:))^0.5);
   u_limit = max(gmmObject.mu(:,j))+3*(max(gmmObject.Sigma(j,j,:))^0.5);
   xmesh_inverse_space{j} = (l_limit:(u_limit-l_limit)/(bin_num-1):u_limit);
end


[~,pdensity{i},xmesh{i}]=kde(currentVar,numMeshPoints);
pdensity{i}(pdensity{i}<0) = 0;

cdensity{i} = cumsum(pdensity{i});
cdensity{i} = (cdensity{i}-min(cdensity{i}))/(max(cdensity{i})-min(cdensity{i})); % scaling the cdensity value to be between [0 1]
    end

[xgrid,ygrid] = meshgrid(xmesh{1}(2:end-1),xmesh{2}(2:end-1));

for k = 1:d
    marginalLogLikelihood_grid{k} = log(pdensity{k}(2:end-1)+eps);
    marginalCDFValues_grid{k} = cdensity{k}(2:end-1);
end
[marg1,marg2] = meshgrid(marginalLogLikelihood_grid{1},marginalLogLikelihood_grid{2});

[xg,yg] = meshgrid(marginalCDFValues_grid{1},marginalCDFValues_grid{2});
inputMatrix = [reshape(xg,numel(xg),1) reshape(yg,numel(yg),1)];


copulaLogLikelihoodVals = gmmCopulaPDF(inputMatrix,gmmObject,xmesh_inverse_space);
Z = reshape(copulaLogLikelihoodVals,size(marg1,1),size(marg1,2));
Z = Z+marg1+marg2;

Z = exp(Z);

plot(givenData(:,1),givenData(:,2),'b.','MarkerSize',3);hold
contour(xgrid,ygrid,Z,40,'EdgeColor',[1 0 0]);

axis tight;

修改方案（实现仅保留轮廓内散点） #

核心思路是提取最外层轮廓的边界坐标，然后判断每个散点是否在轮廓内部，最后只绘制符合条件的点。修改后的完整函数代码如下：

function [xgrid,ygrid,Z,filteredData] = biVariateContourPlotsGMMCopula(givenData,gmmObject,~,numMeshPoints,x_dim,y_dim)

d = 2;
if nargin < 5
    x_dim = 1;
    y_dim = 2;
end

if x_dim == y_dim
    hist(givenData(:,x_dim),10);
    return;
end

numMeshPoints = min(numMeshPoints,256);

% 保留原始数据副本，避免后续修改影响点的判断
originalData = givenData(:,[x_dim y_dim]);
givenData = originalData;
alpha = gmmObject.alpha;
mu = gmmObject.mu(:,[x_dim y_dim]);
sigma = gmmObject.sigma([x_dim y_dim],[x_dim y_dim],:) + 0.005*repmat(eye(d),[1 1 numel(alpha)]);


gmmObject = gmdistribution(mu,sigma,alpha);

bin_num = 256;
for j = 1:2
   l_limit = min(gmmObject.mu(:,j))-3*(max(gmmObject.Sigma(j,j,:))^0.5);
   u_limit = max(gmmObject.mu(:,j))+3*(max(gmmObject.Sigma(j,j,:))^0.5);
   xmesh_inverse_space{j} = (l_limit:(u_limit-l_limit)/(bin_num-1):u_limit);
end

% 补充原代码中缺失的变量遍历逻辑（原代码未定义currentVar）
currentVarList = {givenData(:,1), givenData(:,2)};
for i = 1:d
    [~,pdensity{i},xmesh{i}]=kde(currentVarList{i},numMeshPoints);
    pdensity{i}(pdensity{i}<0) = 0;

    cdensity{i} = cumsum(pdensity{i});
    cdensity{i} = (cdensity{i}-min(cdensity{i}))/(max(cdensity{i})-min(cdensity{i})); % 缩放至[0,1]区间
end

[xgrid,ygrid] = meshgrid(xmesh{1}(2:end-1),xmesh{2}(2:end-1));

for k = 1:d
    marginalLogLikelihood_grid{k} = log(pdensity{k}(2:end-1)+eps);
    marginalCDFValues_grid{k} = cdensity{k}(2:end-1);
end
[marg1,marg2] = meshgrid(marginalLogLikelihood_grid{1},marginalLogLikelihood_grid{2});

[xg,yg] = meshgrid(marginalCDFValues_grid{1},marginalCDFValues_grid{2});
inputMatrix = [reshape(xg,numel(xg),1) reshape(yg,numel(yg),1)];


copulaLogLikelihoodVals = gmmCopulaPDF(inputMatrix,gmmObject,xmesh_inverse_space);
Z = reshape(copulaLogLikelihoodVals,size(marg1,1),size(marg1,2));
Z = Z+marg1+marg2;

Z = exp(Z);

% 获取轮廓数据，提取最外层轮廓（contourc返回的第一个块即为最外层）
contourData = contourc(xgrid,ygrid,Z,40);
idx = 1;
numPoints = contourData(2,idx);
outerContourX = contourData(1,idx+1:idx+numPoints);
outerContourY = contourData(2,idx+1:idx+numPoints);

% 判断每个散点是否在最外层轮廓内部
inFlag = inpolygon(originalData(:,1), originalData(:,2), outerContourX, outerContourY);
% 筛选出轮廓内的点
filteredData = originalData(inFlag,:);

% 绘制筛选后的散点和轮廓
plot(filteredData(:,1),filteredData(:,2),'b.','MarkerSize',3);hold on
contour(xgrid,ygrid,Z,40,'EdgeColor',[1 0 0]);

axis tight;
hold off;

关键修改说明 #

• 补充了原代码中缺失的currentVar遍历逻辑，确保密度计算正常运行。
• 使用contourc提取轮廓数据，直接获取最外层轮廓的坐标集合。
• 通过inpolygon函数逐个判断散点是否在轮廓内部，筛选出符合条件的点后再绘制。
• 新增返回变量filteredData，方便后续直接使用筛选后的数据集。

内容的提问来源于stack exchange，提问作者Arjun Krishnan