[R] Writing an R fn to produce random sequences of binaries XXXX

[Duncan Murdoch]

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 at 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.
>