如何修复斜梯形分布采样输出样本大小

Question

我正在尝试使用逆变换采样来生成倾斜的梯形分布。 输入是斜坡开始和结束的值（a、b、c、d）和样本大小。

a=-3;b=-1;c=1;d=8; 
SampleSize=10e4;
h=2/(d+c-a-b);

然后我计算坡道和平坦部分的长度之比，以获得每个的样本大小：

firstramp=round(((b-a)/(d-a)),3);
flat=round((c-b)/(d-a),3);
secondramp=round((d-c)/(d-a),3);
n1=firstramp*SampleSize; %sample size for first ramp
n3=secondramp*SampleSize; %sample size for second ramp
n2=flat*SampleSize;

最后我从以下代码中得到直方图：

quartile1=h/2*(b-a);
quartile2=1-h/2*(d-c);

y1=linspace(0,quartile1,n1);
y2=linspace(quartile1,quartile2,n2);
y3=linspace(quartile2,1,n3);

%inverse cumulative distribution functions
invcdf1=a+sqrt(2*(b-a)/h)*sqrt(y1);
invcdf2=(a+b)/2+y2/h;
invcdf3=d-sqrt(2*(d-c)/h)*sqrt(1-y3);

distr=[invcdf1 invcdf2 invcdf3];

histogram(distr,100)

但是斜坡和平面分量的采样不相等，看起来像这样：

我通过反复试验解决了这个问题，将坡道的样本量减少了一半：

n1=0.5*firstramp*SampleSize; %sample size for first ramp
n3=0.5*secondramp*SampleSize; %sample size for second ramp
n2=flat*SampleSize;

这使得分布看起来像这样：

但是，这会使输出样本小于输入中给出的样本。

我还尝试了改变坡道和平坦的样本大小的不同组合。 这也有效：

n1=0.75*firstramp*SampleSize; %sample size for first ramp
n3=0.75*secondramp*SampleSize; %sample size for second ramp
n2=1.5*flat*SampleSize;

它增加了输出样本，但仍然没有接近。

任何帮助将不胜感激。

完整代码：

a=-3;b=-1;c=1;d=8; 
SampleSize=10e4;%*1.33333333333333;
h=2/(d+c-a-b);
firstramp=round(((b-a)/(d-a)),3);
flat=round((c-b)/(d-a),3);
secondramp=round((d-c)/(d-a),3);

n1=firstramp*SampleSize; %sample size for first ramp
n3=secondramp*SampleSize; %sample size for second ramp
n2=flat*SampleSize;

quartile1=h/2*(b-a);
quartile2=1-h/2*(d-c);

y1=linspace(0,quartile1,.75*n1);
y2=linspace(quartile1,quartile2,1.5*n2);
y3=linspace(quartile2,1,.75*n3);

%inverse cumulative distribution functions
invcdf1=a+sqrt(2*(b-a)/h)*sqrt(y1);
invcdf2=(a+b)/2+y2/h;
invcdf3=d-sqrt(2*(d-c)/h)*sqrt(1-y3);

distr=[invcdf1 invcdf2 invcdf3];

histogram(distr,100)
%end

Answer 1

我不知道 Matlab，所以我希望其他人会加入这个，但因为没有人在这里做。

如果我正确阅读了您的代码，您所做的并不是倒置。反转是 1-1，即一个统一的输入产生一个结果。您似乎正在使用一种称为“组合方法”的技术。在组合中，整体分布由组件组成，每个组件都可以直接生成。您可以根据它们相对于整体的比例/概率来选择从哪个组件生成。对于密度函数，概率是密度曲线下的面积，因此您的第一个错误是相对于每个分量的宽度而不是使用它们的面积来采样分量。对于您指定的 firstramp、flat 和 secondramp 组件，正确的采样比例分别为 2/13、4/13 和 7/13。第二个错误（相对较小）是将准确的样本大小分配给每个组件。概率为 2/13 并不意味着您的样本中的 2*SampleSize/13 将来自 firstramp，这意味着这是该组件的预期样本大小。随机变量的期望值不一定（甚至可能）是您实际得到的结果。

在伪代码中，组合方法是

generate U ~ Uniform(0,1)
if U <= 2/13:
   generate and return a value from firstramp
else if U <= 6/13:
   generate and return a value from flat
else:
   generate and return a value from secondramp

请注意，由于每个generate 选项都将使用一个或多个统一，并且在选项之间进行选择需要统一U，因此这不是倒置。

如果你想要一个实际的反演，你需要量化你的密度，整合它以获得累积分布函数，然后通过设置F(X) = U 并求解X 来应用反演技术。由于您的分布由不同的分量组成，因此密度和累积密度都是分段函数。

根据两个三角形的面积和平面部分的面积必须为1的要求推导出高度后，我为您的密度得出以下结论：

       | (x + 3) / 13       -3 <= x <= -1
       |
f(x) = | 2 / 13             -1 <= x <= 1
       |
       | 2 * (8 - x) / 91    1 <= x <= 8

整合这个和收集术语产生 CDF：

       | (x + 3)**2 / 26                    -3 <= x <= -1
       |
F(x) = | (2 + x) * 2 / 13                   -1 <= x <= 1
       |
       | 6 / 13 + [49 - (x - 8)**2] / 91     1 <= x <= 8

最后，确定段之间断点处的 F(x) 值并应用反演得到以下伪代码算法：

generate U ~ Uniform(0,1)
if U <= 2 / 13:
    return 2 * sqrt( (13 * U) / 2 ) - 3
else if U <= 6 / 13:
    return (13 * U) / 2 - 2:
else:
    return 8 - sqrt( 91 * (1 - U) )

请注意，这是一个真正的反转。结果是通过生成单个 U 来确定的，并根据它所在的范围以不同的方式对其进行变换。