pearsrnd and pearspdf are not coherent for the Pearson IV distribution?
We are having some trouble using the Pearson 4 implementation of the functions pearsrnd, pearspdf and pearscdf.
We think that pearsrnd does not generate samples from the distribution defined in pearspdf. In particular, we’ve found that pearsrnd over generates observations in the tails of the distribution.
In order to see this, we’ve implemented two simple tests:
1. We’ve compared the PDF of the Pearson 4 distribution (using pearspdf) to the PDF obtained from a simulation using pearsrnd, but we’ve found that the two are very different (even with a very long sample). This problem doesn’t happen for other distributions. In particular, we’ve implemented the same test with the Pearson 6, and the results are satisfactory.
2. In addition, to test if the generatated sample is compatible with its CDF, we’ve applied the function pearscdf to the generated sample. We should have got a uniform distribution, because if X is a random variable and F its CDF, then F(X) is uniform in [0,1]. This is indeed what we’ve got with the Pearson 6, as expected, but not what we’ve got using the Pearson 4.
We’ve conducted a similar experiment using Mathematica, and we get the correct result.
We’ve attached the code we’ve used.
Any idea why we’re getting this result?
Thank you very much.
%% Seed
rng(123)
%% Parameters
%Type IV
mu4 = 5;
sigma4 = 1;
skew4 = 1;
kurt4 = 10;
% Check the distribution type
[ign, type4] = pearsrnd(mu4, sigma4, skew4, kurt4, 1);
disp("This Pearson is of type " + type4)
% Type VI
mu6 = 2;
sigma6 = 1;
skew6 = 2;
kurt6 = 10;
% Check the distribution type
[ign, type6] = pearsrnd(mu6, sigma6, skew6, kurt6, 1);
disp("This Pearson is of type " + type6)
%% Simulate n observations according to both Pearson distributions defined above
n = 100000;
samples4 = pearsrnd(mu4, sigma4, skew4, kurt4, 1, n);
samples6 = pearsrnd(mu6, sigma6, skew6, kurt6, 1, n);
%% First test: we compare the simulated pdf (i.e. histogram of the sample) to the closed form PDF
pdf4 = @(x) pearspdf(x, mu4, sigma4, skew4, kurt4);
pdf6 = @(x) pearspdf(x, mu6, sigma6, skew6, kurt6);
x_ = linspace(-10, 10, 10000);
y4 = pdf4(x_);
y6 = pdf6(x_);
figure(‘Position’, [100, 100, 1000, 400]); % [left, bottom, width, height]
subplot(1, 2, 1);
plot(x_, y4, Color=’red’, LineWidth=2)
title(‘Simulated and closed form PDF for Pearson 4’)
hold on
histogram(samples4, 1000, ‘Normalization’, ‘pdf’, ‘FaceColor’, ‘blue’)
hold off
subplot(1, 2, 2);
plot(x_, y6, Color=’red’, LineWidth=2)
title(‘Simulated and closed form PDF for Pearson 6’)
hold on
histogram(samples6, 1000, ‘Normalization’, ‘pdf’, ‘FaceColor’, ‘blue’)
hold off
sgtitle(‘Test 1: Comparison between Pearson 4 and 6’);
%% Second test: Let X be a random variable and F its CDF, then F(X) is uniform in [0,1].
% We apply the Pearson CDF to the simulated observations and plot its
% histogram. We should obtain something resembling a uniform distribution.
figure(‘Position’, [100, 100, 1000, 400]); % [left, bottom, width, height]
subplot(1, 2, 1);
u4 = pearscdf(samples4, mu4, sigma4, skew4, kurt4);
histogram(u4)
title(‘F(X) for Pearson 4’)
subplot(1, 2, 2);
u6 = pearscdf(samples6, mu6, sigma6, skew6, kurt6);
histogram(u6)
title(‘F(X) for Pearson 6’)
sgtitle(‘Test 2: Comparison between Pearson 4 and 6’);
%% Conclusion
%We get the expected results using the Pearson 6, but we cannot say the
%same for the Pearson 4. Why?We are having some trouble using the Pearson 4 implementation of the functions pearsrnd, pearspdf and pearscdf.
We think that pearsrnd does not generate samples from the distribution defined in pearspdf. In particular, we’ve found that pearsrnd over generates observations in the tails of the distribution.
In order to see this, we’ve implemented two simple tests:
1. We’ve compared the PDF of the Pearson 4 distribution (using pearspdf) to the PDF obtained from a simulation using pearsrnd, but we’ve found that the two are very different (even with a very long sample). This problem doesn’t happen for other distributions. In particular, we’ve implemented the same test with the Pearson 6, and the results are satisfactory.
2. In addition, to test if the generatated sample is compatible with its CDF, we’ve applied the function pearscdf to the generated sample. We should have got a uniform distribution, because if X is a random variable and F its CDF, then F(X) is uniform in [0,1]. This is indeed what we’ve got with the Pearson 6, as expected, but not what we’ve got using the Pearson 4.
We’ve conducted a similar experiment using Mathematica, and we get the correct result.
We’ve attached the code we’ve used.
Any idea why we’re getting this result?
Thank you very much.
%% Seed
rng(123)
%% Parameters
%Type IV
mu4 = 5;
sigma4 = 1;
skew4 = 1;
kurt4 = 10;
% Check the distribution type
[ign, type4] = pearsrnd(mu4, sigma4, skew4, kurt4, 1);
disp("This Pearson is of type " + type4)
% Type VI
mu6 = 2;
sigma6 = 1;
skew6 = 2;
kurt6 = 10;
% Check the distribution type
[ign, type6] = pearsrnd(mu6, sigma6, skew6, kurt6, 1);
disp("This Pearson is of type " + type6)
%% Simulate n observations according to both Pearson distributions defined above
n = 100000;
samples4 = pearsrnd(mu4, sigma4, skew4, kurt4, 1, n);
samples6 = pearsrnd(mu6, sigma6, skew6, kurt6, 1, n);
%% First test: we compare the simulated pdf (i.e. histogram of the sample) to the closed form PDF
pdf4 = @(x) pearspdf(x, mu4, sigma4, skew4, kurt4);
pdf6 = @(x) pearspdf(x, mu6, sigma6, skew6, kurt6);
x_ = linspace(-10, 10, 10000);
y4 = pdf4(x_);
y6 = pdf6(x_);
figure(‘Position’, [100, 100, 1000, 400]); % [left, bottom, width, height]
subplot(1, 2, 1);
plot(x_, y4, Color=’red’, LineWidth=2)
title(‘Simulated and closed form PDF for Pearson 4’)
hold on
histogram(samples4, 1000, ‘Normalization’, ‘pdf’, ‘FaceColor’, ‘blue’)
hold off
subplot(1, 2, 2);
plot(x_, y6, Color=’red’, LineWidth=2)
title(‘Simulated and closed form PDF for Pearson 6’)
hold on
histogram(samples6, 1000, ‘Normalization’, ‘pdf’, ‘FaceColor’, ‘blue’)
hold off
sgtitle(‘Test 1: Comparison between Pearson 4 and 6’);
%% Second test: Let X be a random variable and F its CDF, then F(X) is uniform in [0,1].
% We apply the Pearson CDF to the simulated observations and plot its
% histogram. We should obtain something resembling a uniform distribution.
figure(‘Position’, [100, 100, 1000, 400]); % [left, bottom, width, height]
subplot(1, 2, 1);
u4 = pearscdf(samples4, mu4, sigma4, skew4, kurt4);
histogram(u4)
title(‘F(X) for Pearson 4’)
subplot(1, 2, 2);
u6 = pearscdf(samples6, mu6, sigma6, skew6, kurt6);
histogram(u6)
title(‘F(X) for Pearson 6’)
sgtitle(‘Test 2: Comparison between Pearson 4 and 6’);
%% Conclusion
%We get the expected results using the Pearson 6, but we cannot say the
%same for the Pearson 4. Why? We are having some trouble using the Pearson 4 implementation of the functions pearsrnd, pearspdf and pearscdf.
We think that pearsrnd does not generate samples from the distribution defined in pearspdf. In particular, we’ve found that pearsrnd over generates observations in the tails of the distribution.
In order to see this, we’ve implemented two simple tests:
1. We’ve compared the PDF of the Pearson 4 distribution (using pearspdf) to the PDF obtained from a simulation using pearsrnd, but we’ve found that the two are very different (even with a very long sample). This problem doesn’t happen for other distributions. In particular, we’ve implemented the same test with the Pearson 6, and the results are satisfactory.
2. In addition, to test if the generatated sample is compatible with its CDF, we’ve applied the function pearscdf to the generated sample. We should have got a uniform distribution, because if X is a random variable and F its CDF, then F(X) is uniform in [0,1]. This is indeed what we’ve got with the Pearson 6, as expected, but not what we’ve got using the Pearson 4.
We’ve conducted a similar experiment using Mathematica, and we get the correct result.
We’ve attached the code we’ve used.
Any idea why we’re getting this result?
Thank you very much.
%% Seed
rng(123)
%% Parameters
%Type IV
mu4 = 5;
sigma4 = 1;
skew4 = 1;
kurt4 = 10;
% Check the distribution type
[ign, type4] = pearsrnd(mu4, sigma4, skew4, kurt4, 1);
disp("This Pearson is of type " + type4)
% Type VI
mu6 = 2;
sigma6 = 1;
skew6 = 2;
kurt6 = 10;
% Check the distribution type
[ign, type6] = pearsrnd(mu6, sigma6, skew6, kurt6, 1);
disp("This Pearson is of type " + type6)
%% Simulate n observations according to both Pearson distributions defined above
n = 100000;
samples4 = pearsrnd(mu4, sigma4, skew4, kurt4, 1, n);
samples6 = pearsrnd(mu6, sigma6, skew6, kurt6, 1, n);
%% First test: we compare the simulated pdf (i.e. histogram of the sample) to the closed form PDF
pdf4 = @(x) pearspdf(x, mu4, sigma4, skew4, kurt4);
pdf6 = @(x) pearspdf(x, mu6, sigma6, skew6, kurt6);
x_ = linspace(-10, 10, 10000);
y4 = pdf4(x_);
y6 = pdf6(x_);
figure(‘Position’, [100, 100, 1000, 400]); % [left, bottom, width, height]
subplot(1, 2, 1);
plot(x_, y4, Color=’red’, LineWidth=2)
title(‘Simulated and closed form PDF for Pearson 4’)
hold on
histogram(samples4, 1000, ‘Normalization’, ‘pdf’, ‘FaceColor’, ‘blue’)
hold off
subplot(1, 2, 2);
plot(x_, y6, Color=’red’, LineWidth=2)
title(‘Simulated and closed form PDF for Pearson 6’)
hold on
histogram(samples6, 1000, ‘Normalization’, ‘pdf’, ‘FaceColor’, ‘blue’)
hold off
sgtitle(‘Test 1: Comparison between Pearson 4 and 6’);
%% Second test: Let X be a random variable and F its CDF, then F(X) is uniform in [0,1].
% We apply the Pearson CDF to the simulated observations and plot its
% histogram. We should obtain something resembling a uniform distribution.
figure(‘Position’, [100, 100, 1000, 400]); % [left, bottom, width, height]
subplot(1, 2, 1);
u4 = pearscdf(samples4, mu4, sigma4, skew4, kurt4);
histogram(u4)
title(‘F(X) for Pearson 4’)
subplot(1, 2, 2);
u6 = pearscdf(samples6, mu6, sigma6, skew6, kurt6);
histogram(u6)
title(‘F(X) for Pearson 6’)
sgtitle(‘Test 2: Comparison between Pearson 4 and 6’);
%% Conclusion
%We get the expected results using the Pearson 6, but we cannot say the
%same for the Pearson 4. Why? pearson, distribution, statistics, random number generator, type iv MATLAB Answers — New Questions
