On 13-10-15 7:55 AM, Dan Abner wrote:
Hi all,
I am attempting to write an R fn to produce a random sequence of 8 binaries
(4 ordered pairs). I would like to parmeterize the fn so that it takes 2
arguments: 1 that is an overall probability of a 1 vs. 0 and the other
which controls the likelihood of 1s on the 1st vs. the 2nd element of the 4
ordered pairs.
Here are some examples:
*Voiced Probability Parm* *1st Element Probability Parm* *Average N of 1s
* *Average N 1st Element 1s* * * *1* *1.5* *2* *2.5* *3* *3.5* *4*
*4.5*
0.5 1 4 4 1 0 1 0 1 0 1 0 0.5 0 4 0 0 1 0 1 0 1 0 1 0.5 0.75 4 3 1 0 0 1 1
0 1 0 0.5 0.25 4 1 0 1 0 1 1 0 0 1 1 1 8 4 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0
0 0 0 0 0 0.5 0.5 4 2 1 0 0 1 1 0 0 1
Can anyone suggest an algorithm for this? I have a 1st draft of the fn
definition, but my algorithm is not correct as I am not getting the average
number of 1s (averaged over 10,000 sequences) correct.
A simple algorithm is to sample the 1st bit using rbinom(), then sample
the second bit, again using rbinom, but with probability conditional on
the value of the first bit. But that requires a joint distribution, and
you haven't given us enough information to construct one.
For example, here's a model where bit 2 is either independent of bit 1
or forced to be equal to it, with theta being the probability of being
forced to be equal. This probably isn't the model you want, but you
would just modify the bit2 line to give the correct conditional
distribution based on what you do want.
p <- 0.7 # marginal prob of a 1
theta <- 0.5 # ratio of times bit 2 is equal to bit 1, versus
# independent
bit1 <- rbinom(4, size=1, p)
bit2 <- rbinom(4, size=1, theta*bit1 + (1-theta)*p)
cbind(bit1, bit2)
Duncan Murdoch
Thanks,
Dan
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
