Pymc3 model where the results of a switch are directly observed
I have just started learning pymc3 so I might be thinking about this completely the wrong way.
Assume that we observe a vector of 10 booleans.
The process of interest generates (observed) booleans with a Bernoulli distribution with a parameter theta1. So I define a Beta prior over theta1 and define a variable with length 10 that is a sample from Bernoulli(theta1).
However, this true sample is disturbed by sometimes switching the true data to 0, with a probability theta2. So I define a switch to 0 with a probability Bernoulli(theta2).
The switched values are the observed ones. I am not sure how to tell the model that I observed the switched variables, i.e. I am not sure how to fit the model to the observed data.
This is what I have for now, and I am kind of stuck:
# observed data (already switched)
observed_data = np.random.binomial(1, 0.5, size=10)
with pm.Model() as skeptic_model:
# uniform probability of the bernoulli parameter
true_model_prior = pm.Beta("true_model_prior", 1, 1)
true_data = pm.Bernoulli("true_data", p=true_model_prior, shape=data.shape)
disturbed_data = pm.math.switch(pm.Bernoulli("disturbed", 0.1), true_data, 0)
python-3.x pymc3
asked Nov 10 at 17:07
whatamess
639
add a comment |
I have just started learning pymc3 so I might be thinking about this completely the wrong way.
Assume that we observe a vector of 10 booleans.
The process of interest generates (observed) booleans with a Bernoulli distribution with a parameter theta1. So I define a Beta prior over theta1 and define a variable with length 10 that is a sample from Bernoulli(theta1).
However, this true sample is disturbed by sometimes switching the true data to 0, with a probability theta2. So I define a switch to 0 with a probability Bernoulli(theta2).
The switched values are the observed ones. I am not sure how to tell the model that I observed the switched variables, i.e. I am not sure how to fit the model to the observed data.
This is what I have for now, and I am kind of stuck:
# observed data (already switched)
observed_data = np.random.binomial(1, 0.5, size=10)
with pm.Model() as skeptic_model:
# uniform probability of the bernoulli parameter
true_model_prior = pm.Beta("true_model_prior", 1, 1)
true_data = pm.Bernoulli("true_data", p=true_model_prior, shape=data.shape)
disturbed_data = pm.math.switch(pm.Bernoulli("disturbed", 0.1), true_data, 0)
python-3.x pymc3
asked Nov 10 at 17:07
whatamess
639
Before usingpymc3, you have to create that second array, the one with the switched values. How do you do it?
– David
Nov 13 at 6:41
add a comment |
I have just started learning pymc3 so I might be thinking about this completely the wrong way.
Assume that we observe a vector of 10 booleans.
The process of interest generates (observed) booleans with a Bernoulli distribution with a parameter theta1. So I define a Beta prior over theta1 and define a variable with length 10 that is a sample from Bernoulli(theta1).
However, this true sample is disturbed by sometimes switching the true data to 0, with a probability theta2. So I define a switch to 0 with a probability Bernoulli(theta2).
The switched values are the observed ones. I am not sure how to tell the model that I observed the switched variables, i.e. I am not sure how to fit the model to the observed data.
This is what I have for now, and I am kind of stuck:
# observed data (already switched)
observed_data = np.random.binomial(1, 0.5, size=10)
with pm.Model() as skeptic_model:
# uniform probability of the bernoulli parameter
true_model_prior = pm.Beta("true_model_prior", 1, 1)
true_data = pm.Bernoulli("true_data", p=true_model_prior, shape=data.shape)
disturbed_data = pm.math.switch(pm.Bernoulli("disturbed", 0.1), true_data, 0)
python-3.x pymc3
asked Nov 10 at 17:07
whatamess
639
I have just started learning pymc3 so I might be thinking about this completely the wrong way.
Assume that we observe a vector of 10 booleans.
The process of interest generates (observed) booleans with a Bernoulli distribution with a parameter theta1. So I define a Beta prior over theta1 and define a variable with length 10 that is a sample from Bernoulli(theta1).
However, this true sample is disturbed by sometimes switching the true data to 0, with a probability theta2. So I define a switch to 0 with a probability Bernoulli(theta2).
The switched values are the observed ones. I am not sure how to tell the model that I observed the switched variables, i.e. I am not sure how to fit the model to the observed data.
This is what I have for now, and I am kind of stuck:
# observed data (already switched)
observed_data = np.random.binomial(1, 0.5, size=10)
with pm.Model() as skeptic_model:
# uniform probability of the bernoulli parameter
true_model_prior = pm.Beta("true_model_prior", 1, 1)
true_data = pm.Bernoulli("true_data", p=true_model_prior, shape=data.shape)
disturbed_data = pm.math.switch(pm.Bernoulli("disturbed", 0.1), true_data, 0)
python-3.x pymc3
python-3.x pymc3
asked Nov 10 at 17:07
whatamess
639
asked Nov 10 at 17:07
whatamess
639
asked Nov 10 at 17:07
whatamess
639
asked Nov 10 at 17:07
whatamess
639
asked Nov 10 at 17:07
whatamess
639
639
Before usingpymc3, you have to create that second array, the one with the switched values. How do you do it?
– David
Nov 13 at 6:41
add a comment |
Before usingpymc3, you have to create that second array, the one with the switched values. How do you do it?
– David
Nov 13 at 6:41
Before using
pymc3, you have to create that second array, the one with the switched values. How do you do it?– David
Nov 13 at 6:41
Before using
pymc3, you have to create that second array, the one with the switched values. How do you do it?– David
Nov 13 at 6:41
add a comment |
1 Answer
1
active
oldest
votes
Your model can be reframed as a product of Bernoulli random variables, and therefore as a single Bernoulli random variable with a multiplicative p. Namely, the following model is equivalent to yours:
# observed data (already considered zero-inflated)
Y = np.random.binomial(1, 0.5, size=10)
with pm.Model() as zero_inflated_beta_bernoulli:
# true_model_prior
p = pm.Beta('p', alpha=1, beta=1)
# dropout rate
d = 0.1
# disturbed_data;
y = pm.Bernoulli('y', p = (1-d)*p, observed=Y)
You could let the dropout rate also be a random variable,
# dropout rate
d = pm.Beta('d', mu=0.1, sd=0.02)
However, it should be noted that this model really can't distinguish between dropouts and original outcomes, so the posteriors are sensitive to the priors.
answered Nov 13 at 18:31
merv
24.7k671109
add a comment |
Your Answer
StackExchange.ifUsing("editor", function () {
StackExchange.using("externalEditor", function () {
StackExchange.using("snippets", function () {
StackExchange.snippets.init();
});
});
}, "code-snippets");
StackExchange.ready(function() {
var channelOptions = {
tags: "".split(" "),
id: "1"
};
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
},
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%2fstackoverflow.com%2fquestions%2f53241342%2fpymc3-model-where-the-results-of-a-switch-are-directly-observed%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
Your model can be reframed as a product of Bernoulli random variables, and therefore as a single Bernoulli random variable with a multiplicative p. Namely, the following model is equivalent to yours:
# observed data (already considered zero-inflated)
Y = np.random.binomial(1, 0.5, size=10)
with pm.Model() as zero_inflated_beta_bernoulli:
# true_model_prior
p = pm.Beta('p', alpha=1, beta=1)
# dropout rate
d = 0.1
# disturbed_data;
y = pm.Bernoulli('y', p = (1-d)*p, observed=Y)
You could let the dropout rate also be a random variable,
# dropout rate
d = pm.Beta('d', mu=0.1, sd=0.02)
However, it should be noted that this model really can't distinguish between dropouts and original outcomes, so the posteriors are sensitive to the priors.
answered Nov 13 at 18:31
merv
24.7k671109
add a comment |
Your model can be reframed as a product of Bernoulli random variables, and therefore as a single Bernoulli random variable with a multiplicative p. Namely, the following model is equivalent to yours:
# observed data (already considered zero-inflated)
Y = np.random.binomial(1, 0.5, size=10)
with pm.Model() as zero_inflated_beta_bernoulli:
# true_model_prior
p = pm.Beta('p', alpha=1, beta=1)
# dropout rate
d = 0.1
# disturbed_data;
y = pm.Bernoulli('y', p = (1-d)*p, observed=Y)
You could let the dropout rate also be a random variable,
# dropout rate
d = pm.Beta('d', mu=0.1, sd=0.02)
However, it should be noted that this model really can't distinguish between dropouts and original outcomes, so the posteriors are sensitive to the priors.
answered Nov 13 at 18:31
merv
24.7k671109
add a comment |
Your model can be reframed as a product of Bernoulli random variables, and therefore as a single Bernoulli random variable with a multiplicative p. Namely, the following model is equivalent to yours:
# observed data (already considered zero-inflated)
Y = np.random.binomial(1, 0.5, size=10)
with pm.Model() as zero_inflated_beta_bernoulli:
# true_model_prior
p = pm.Beta('p', alpha=1, beta=1)
# dropout rate
d = 0.1
# disturbed_data;
y = pm.Bernoulli('y', p = (1-d)*p, observed=Y)
You could let the dropout rate also be a random variable,
# dropout rate
d = pm.Beta('d', mu=0.1, sd=0.02)
However, it should be noted that this model really can't distinguish between dropouts and original outcomes, so the posteriors are sensitive to the priors.
answered Nov 13 at 18:31
merv
24.7k671109
Your model can be reframed as a product of Bernoulli random variables, and therefore as a single Bernoulli random variable with a multiplicative p. Namely, the following model is equivalent to yours:
# observed data (already considered zero-inflated)
Y = np.random.binomial(1, 0.5, size=10)
with pm.Model() as zero_inflated_beta_bernoulli:
# true_model_prior
p = pm.Beta('p', alpha=1, beta=1)
# dropout rate
d = 0.1
# disturbed_data;
y = pm.Bernoulli('y', p = (1-d)*p, observed=Y)
You could let the dropout rate also be a random variable,
# dropout rate
d = pm.Beta('d', mu=0.1, sd=0.02)
However, it should be noted that this model really can't distinguish between dropouts and original outcomes, so the posteriors are sensitive to the priors.
answered Nov 13 at 18:31
merv
24.7k671109
edited Nov 14 at 22:59
answered Nov 13 at 18:31
merv
24.7k671109
answered Nov 13 at 18:31
merv
24.7k671109
answered Nov 13 at 18:31
merv
24.7k671109
24.7k671109
add a comment |
add a comment |
Thanks for contributing an answer to Stack Overflow!
- 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.
To learn more, see our tips on writing great answers.
Some of your past answers have not been well-received, and you're in danger of being blocked from answering.
Please pay close attention to the following guidance:
- 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.
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%2fstackoverflow.com%2fquestions%2f53241342%2fpymc3-model-where-the-results-of-a-switch-are-directly-observed%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
Before using
pymc3, you have to create that second array, the one with the switched values. How do you do it?– David
Nov 13 at 6:41