simulate random variable from pdf (discrete case)
$begingroup$
I am trying to simulate n random discrete variable which has the following
pmf
$P(X = k) = (1-p)^2kp^{k-1}$
I am thinking about using the inverse transform sampling method and I am trying to find the cdf.
$P(X le k) = 1 - P(X > k) = 1- P(X ge k+1) = sum_{x=k+1}^{infty} (1-p)^2xp^{x-1} = 1-[ (1-p)^2sum_{x=k+1}^{infty}xp^{x-1}]$
= $1-(1-p)^2[(k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...] $
What I have done so far :
Let $S = (k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...$
$pS =$ $(k+1)p^{k+1} + (k+2)p^{k+2}+ (k+3)p^{k+3} + ...$
then $(1-p)S = kp^k +p^k + p^{k+1} + p^{k+2}+p^{k+3}+...$
(geometric series )
$(1-p)S = kp^k + frac{p^k}{1-p}$
$ F(x) = 1- (xp^x(1-p) + p^x)$
Then I'm trying to find the inverse :
$1- (xp^x(1-p) + p^x) < U le 1- ((x+1)p^{x+1}(1-p) + p^{x+1})$
$xp^x(1-p) + p^x < 1-U < (x+1)p^{x+1}(1-p) + p^{x+1}$
I'm stuck here ...
Any help or hint will be appreciated !
calculus probability statistics matlab
asked Dec 5 '18 at 8:14
TatariaTataria
277
$endgroup$
add a comment |
$begingroup$
I am trying to simulate n random discrete variable which has the following
pmf
$P(X = k) = (1-p)^2kp^{k-1}$
I am thinking about using the inverse transform sampling method and I am trying to find the cdf.
$P(X le k) = 1 - P(X > k) = 1- P(X ge k+1) = sum_{x=k+1}^{infty} (1-p)^2xp^{x-1} = 1-[ (1-p)^2sum_{x=k+1}^{infty}xp^{x-1}]$
= $1-(1-p)^2[(k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...] $
What I have done so far :
Let $S = (k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...$
$pS =$ $(k+1)p^{k+1} + (k+2)p^{k+2}+ (k+3)p^{k+3} + ...$
then $(1-p)S = kp^k +p^k + p^{k+1} + p^{k+2}+p^{k+3}+...$
(geometric series )
$(1-p)S = kp^k + frac{p^k}{1-p}$
$ F(x) = 1- (xp^x(1-p) + p^x)$
Then I'm trying to find the inverse :
$1- (xp^x(1-p) + p^x) < U le 1- ((x+1)p^{x+1}(1-p) + p^{x+1})$
$xp^x(1-p) + p^x < 1-U < (x+1)p^{x+1}(1-p) + p^{x+1}$
I'm stuck here ...
Any help or hint will be appreciated !
calculus probability statistics matlab
asked Dec 5 '18 at 8:14
TatariaTataria
277
$endgroup$
1
$begingroup$
I do not think there is a nice inverse for the CDF. For discrete random variable, as the CDF is just a step function, you may simply do a summation of the pmf to obtain that, and use that to generate. For example, in your case the support is ${1, 2, 3, ldots}$, you generate $U sim text{Uniform}(0, 1)$, then check: If $U < f_X(1)$, where $f_X$ is the pmf of $X$, then assign $X$ = 1. Else check if $U < f_X(1) + f_X(2)$, then assign $X = 2$, and so on. The speed is not too bad if $f(x)$ is mainly dominated in the first few terms, and the tail rapidly converge.
$endgroup$
– BGM
Dec 5 '18 at 9:12
$begingroup$
I guess I have no choice
$endgroup$
– Tataria
Dec 5 '18 at 9:19
add a comment |
$begingroup$
I am trying to simulate n random discrete variable which has the following
pmf
$P(X = k) = (1-p)^2kp^{k-1}$
I am thinking about using the inverse transform sampling method and I am trying to find the cdf.
$P(X le k) = 1 - P(X > k) = 1- P(X ge k+1) = sum_{x=k+1}^{infty} (1-p)^2xp^{x-1} = 1-[ (1-p)^2sum_{x=k+1}^{infty}xp^{x-1}]$
= $1-(1-p)^2[(k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...] $
What I have done so far :
Let $S = (k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...$
$pS =$ $(k+1)p^{k+1} + (k+2)p^{k+2}+ (k+3)p^{k+3} + ...$
then $(1-p)S = kp^k +p^k + p^{k+1} + p^{k+2}+p^{k+3}+...$
(geometric series )
$(1-p)S = kp^k + frac{p^k}{1-p}$
$ F(x) = 1- (xp^x(1-p) + p^x)$
Then I'm trying to find the inverse :
$1- (xp^x(1-p) + p^x) < U le 1- ((x+1)p^{x+1}(1-p) + p^{x+1})$
$xp^x(1-p) + p^x < 1-U < (x+1)p^{x+1}(1-p) + p^{x+1}$
I'm stuck here ...
Any help or hint will be appreciated !
calculus probability statistics matlab
asked Dec 5 '18 at 8:14
TatariaTataria
277
$endgroup$
I am trying to simulate n random discrete variable which has the following
pmf
$P(X = k) = (1-p)^2kp^{k-1}$
I am thinking about using the inverse transform sampling method and I am trying to find the cdf.
$P(X le k) = 1 - P(X > k) = 1- P(X ge k+1) = sum_{x=k+1}^{infty} (1-p)^2xp^{x-1} = 1-[ (1-p)^2sum_{x=k+1}^{infty}xp^{x-1}]$
= $1-(1-p)^2[(k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...] $
What I have done so far :
Let $S = (k+1)p^k + (k+2)p^{k+1}+ (k+3)p^{k+2} + ...$
$pS =$ $(k+1)p^{k+1} + (k+2)p^{k+2}+ (k+3)p^{k+3} + ...$
then $(1-p)S = kp^k +p^k + p^{k+1} + p^{k+2}+p^{k+3}+...$
(geometric series )
$(1-p)S = kp^k + frac{p^k}{1-p}$
$ F(x) = 1- (xp^x(1-p) + p^x)$
Then I'm trying to find the inverse :
$1- (xp^x(1-p) + p^x) < U le 1- ((x+1)p^{x+1}(1-p) + p^{x+1})$
$xp^x(1-p) + p^x < 1-U < (x+1)p^{x+1}(1-p) + p^{x+1}$
I'm stuck here ...
Any help or hint will be appreciated !
calculus probability statistics matlab
calculus probability statistics matlab
asked Dec 5 '18 at 8:14
TatariaTataria
277
asked Dec 5 '18 at 8:14
TatariaTataria
277
asked Dec 5 '18 at 8:14
TatariaTataria
277
asked Dec 5 '18 at 8:14
TatariaTataria
277
asked Dec 5 '18 at 8:14
TatariaTataria
277
277
1
$begingroup$
I do not think there is a nice inverse for the CDF. For discrete random variable, as the CDF is just a step function, you may simply do a summation of the pmf to obtain that, and use that to generate. For example, in your case the support is ${1, 2, 3, ldots}$, you generate $U sim text{Uniform}(0, 1)$, then check: If $U < f_X(1)$, where $f_X$ is the pmf of $X$, then assign $X$ = 1. Else check if $U < f_X(1) + f_X(2)$, then assign $X = 2$, and so on. The speed is not too bad if $f(x)$ is mainly dominated in the first few terms, and the tail rapidly converge.
$endgroup$
– BGM
Dec 5 '18 at 9:12
$begingroup$
I guess I have no choice
$endgroup$
– Tataria
Dec 5 '18 at 9:19
add a comment |
1
$begingroup$
I do not think there is a nice inverse for the CDF. For discrete random variable, as the CDF is just a step function, you may simply do a summation of the pmf to obtain that, and use that to generate. For example, in your case the support is ${1, 2, 3, ldots}$, you generate $U sim text{Uniform}(0, 1)$, then check: If $U < f_X(1)$, where $f_X$ is the pmf of $X$, then assign $X$ = 1. Else check if $U < f_X(1) + f_X(2)$, then assign $X = 2$, and so on. The speed is not too bad if $f(x)$ is mainly dominated in the first few terms, and the tail rapidly converge.
$endgroup$
– BGM
Dec 5 '18 at 9:12
$begingroup$
I guess I have no choice
$endgroup$
– Tataria
Dec 5 '18 at 9:19
1
1
$begingroup$
I do not think there is a nice inverse for the CDF. For discrete random variable, as the CDF is just a step function, you may simply do a summation of the pmf to obtain that, and use that to generate. For example, in your case the support is ${1, 2, 3, ldots}$, you generate $U sim text{Uniform}(0, 1)$, then check: If $U < f_X(1)$, where $f_X$ is the pmf of $X$, then assign $X$ = 1. Else check if $U < f_X(1) + f_X(2)$, then assign $X = 2$, and so on. The speed is not too bad if $f(x)$ is mainly dominated in the first few terms, and the tail rapidly converge.
$endgroup$
– BGM
Dec 5 '18 at 9:12
$begingroup$
I do not think there is a nice inverse for the CDF. For discrete random variable, as the CDF is just a step function, you may simply do a summation of the pmf to obtain that, and use that to generate. For example, in your case the support is ${1, 2, 3, ldots}$, you generate $U sim text{Uniform}(0, 1)$, then check: If $U < f_X(1)$, where $f_X$ is the pmf of $X$, then assign $X$ = 1. Else check if $U < f_X(1) + f_X(2)$, then assign $X = 2$, and so on. The speed is not too bad if $f(x)$ is mainly dominated in the first few terms, and the tail rapidly converge.
$endgroup$
– BGM
Dec 5 '18 at 9:12
$begingroup$
I guess I have no choice
$endgroup$
– Tataria
Dec 5 '18 at 9:19
$begingroup$
I guess I have no choice
$endgroup$
– Tataria
Dec 5 '18 at 9:19
add a comment |
1 Answer
1
active
oldest
votes
$begingroup$
There are several ways you can randomly sample from a discrete distribution, unfortunately inverting the CDF is not one of them. The plot below was generated using the Alias Method, it is particularly efficient if you know how to use binary search trees, otherwise a simple implementation of argmin, argmax will work
The blue bars are a histogram of $400$ samples generated with the alias method, the red points are simply
$$
P(X = x) = (1 - p)^2 x p^{x - 1}
$$
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
$endgroup$
add a comment |
Your Answer
StackExchange.ifUsing("editor", function () {
return StackExchange.using("mathjaxEditing", function () {
StackExchange.MarkdownEditor.creationCallbacks.add(function (editor, postfix) {
StackExchange.mathjaxEditing.prepareWmdForMathJax(editor, postfix, [["$", "$"], ["\\(","\\)"]]);
});
});
}, "mathjax-editing");
StackExchange.ready(function() {
var channelOptions = {
tags: "".split(" "),
id: "69"
};
initTagRenderer("".split(" "), "".split(" "), channelOptions);
StackExchange.using("externalEditor", function() {
// Have to fire editor after snippets, if snippets enabled
if (StackExchange.settings.snippets.snippetsEnabled) {
StackExchange.using("snippets", function() {
createEditor();
});
}
else {
createEditor();
}
});
function createEditor() {
StackExchange.prepareEditor({
heartbeatType: 'answer',
autoActivateHeartbeat: false,
convertImagesToLinks: true,
noModals: true,
showLowRepImageUploadWarning: true,
reputationToPostImages: 10,
bindNavPrevention: true,
postfix: "",
imageUploader: {
brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
allowUrls: true
},
noCode: true, onDemand: true,
discardSelector: ".discard-answer"
,immediatelyShowMarkdownHelp:true
});
}
});
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fmath.stackexchange.com%2fquestions%2f3026799%2fsimulate-random-variable-from-pdf-discrete-case%23new-answer', 'question_page');
}
);
Post as a guest
Required, but never shown
1 Answer
1
active
oldest
votes
1 Answer
1
active
oldest
votes
active
oldest
votes
active
oldest
votes
$begingroup$
There are several ways you can randomly sample from a discrete distribution, unfortunately inverting the CDF is not one of them. The plot below was generated using the Alias Method, it is particularly efficient if you know how to use binary search trees, otherwise a simple implementation of argmin, argmax will work
The blue bars are a histogram of $400$ samples generated with the alias method, the red points are simply
$$
P(X = x) = (1 - p)^2 x p^{x - 1}
$$
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
$endgroup$
add a comment |
$begingroup$
There are several ways you can randomly sample from a discrete distribution, unfortunately inverting the CDF is not one of them. The plot below was generated using the Alias Method, it is particularly efficient if you know how to use binary search trees, otherwise a simple implementation of argmin, argmax will work
The blue bars are a histogram of $400$ samples generated with the alias method, the red points are simply
$$
P(X = x) = (1 - p)^2 x p^{x - 1}
$$
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
$endgroup$
add a comment |
$begingroup$
There are several ways you can randomly sample from a discrete distribution, unfortunately inverting the CDF is not one of them. The plot below was generated using the Alias Method, it is particularly efficient if you know how to use binary search trees, otherwise a simple implementation of argmin, argmax will work
The blue bars are a histogram of $400$ samples generated with the alias method, the red points are simply
$$
P(X = x) = (1 - p)^2 x p^{x - 1}
$$
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
$endgroup$
There are several ways you can randomly sample from a discrete distribution, unfortunately inverting the CDF is not one of them. The plot below was generated using the Alias Method, it is particularly efficient if you know how to use binary search trees, otherwise a simple implementation of argmin, argmax will work
The blue bars are a histogram of $400$ samples generated with the alias method, the red points are simply
$$
P(X = x) = (1 - p)^2 x p^{x - 1}
$$
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
answered Dec 5 '18 at 11:09
caveraccaverac
14.2k21130
14.2k21130
add a comment |
add a comment |
Thanks for contributing an answer to Mathematics Stack Exchange!
- Please be sure to answer the question. Provide details and share your research!
But avoid …
- Asking for help, clarification, or responding to other answers.
- Making statements based on opinion; back them up with references or personal experience.
Use MathJax to format equations. MathJax reference.
To learn more, see our tips on writing great answers.
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fmath.stackexchange.com%2fquestions%2f3026799%2fsimulate-random-variable-from-pdf-discrete-case%23new-answer', 'question_page');
}
);
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
1
$begingroup$
I do not think there is a nice inverse for the CDF. For discrete random variable, as the CDF is just a step function, you may simply do a summation of the pmf to obtain that, and use that to generate. For example, in your case the support is ${1, 2, 3, ldots}$, you generate $U sim text{Uniform}(0, 1)$, then check: If $U < f_X(1)$, where $f_X$ is the pmf of $X$, then assign $X$ = 1. Else check if $U < f_X(1) + f_X(2)$, then assign $X = 2$, and so on. The speed is not too bad if $f(x)$ is mainly dominated in the first few terms, and the tail rapidly converge.
$endgroup$
– BGM
Dec 5 '18 at 9:12
$begingroup$
I guess I have no choice
$endgroup$
– Tataria
Dec 5 '18 at 9:19