On Tue, 15 Dec 2009, li li wrote:
Hi:
Thank you for your reply. I realize that I did not state the problem
clearly before.
Here is the problem again.
Let X = (X1, X2, ..., Xn) and Y=(Y1,Y2,...,Yn). Xi's can be 0 or 1.
When Xi=1, Yi is distributed as N(2,1). When Xi=0, Yi is distributed as
N(0,1).
There are 2^n possible X values. For example, x=(0,0,...,0) , i.e. all xi is
0.
The summand is:
beta(a+sum(xi), b+n-sum(xi) ) * [ (1-x1) * dnorm(y1, 0, 1) + x1 * dnorm(y1,
2, 1) ] *
[(1-x2) * dnorm(y2, 0, 1) + x2 * dnorm(y2, 2, 1) ] * ...* [ (1-xn) *
dnorm(yn, 0, 1) + xn * dnorm(yn, 2, 1) ]
where sum(xi)=x1 + x2 + ... + xn.
For example, when x=(0,0,...,0) , i.e. all xi is 0. The summand is
beta(a, b+n) * dnorm(y1, 0, 1) * dnorm(y2, 0, 1) * ...* dnorm(yn, 0, 1)
I want to take the sum of the summand over all possible X values. There are
2^n of them.
I am wondering whether there is a function in R that can take care of this
kind of sum.
A function that does exactly that? I do not know of one, but have not
searched very hard.
If I understood what you want (and "commented, minimal, self-contained,
reproducible code" was requested and would have helped me), a simple one
line command would give the result.
For
set.seed(1)
n <- 10
y <- rnorm(n)
x <- rbinom(10,1,.5)
a <- 1
b <- 2
The answer is given by:
# [ one line command omitted ]
[1] 1.298587e-08
You will need to understand "vectorization" and "recycling rules" to do
this neatly in one line.
It will help you to study the sections of the manuals that pertain to
those topics (and perhaps the mail ist archive) as well as to study the
help and examples for
?dnorm
?beta
?expand.grid
?matrix
?rowSums
?as.matrix
And here is a helpful hint
xmat <- as.matrix( expand.grid( rep( list(0:1), n ) )
will be useful.
If you need more help, be sure to review the Posting Guide and run the
function help.request() to help you better frame your question.
HTH,
Chuck
p.s. Is this homework? If so, ask your instructor for help first.
Thank you
2009/12/15 Dennis Murphy <djmu...@gmail.com>
Hi:
Your expression makes no sense in that I don't see where the summation can
obtain 'over all 2^n possible values...' the way you wrote it. How about
this?
Let x = (x1, x2, ..., xn). The summand can be then be expressed in R as
beta(a + sum(x), b + n - sum(x)) * ((1 - x) * dnorm(x, 0, 1) + x *
dnorm(x, 2, 1)) .
The beta function term is a scalar, the second expression is a vector.
Do you want the sum of the products of the vector elements? If so, it
should be as easy as
beta(a + sum(x), b + n - sum(x)) * sum((1 - x) * dnorm(x, 0, 1) + x *
dnorm(x, 2, 1))
If, on the other hand, you meant to write the summand as
beta(a + x, b + 1 - x) * ((1 - x) * dnorm(x, 0, 1) + x * dnorm(x, 2, 1)) ,
then the sum is
sum(beta(a + x, b + 1 - x) * ((1 - x) * dnorm(x, 0, 1) + x * dnorm(x, 2,
1))
More below...HTH
Dennis
On Mon, Dec 14, 2009 at 6:38 PM, li li <hannah....@gmail.com> wrote:
Hello,
Can anyone give me some suggestion in term of calculating the sum below.
Is there a function in R that can help doing it faster?
x1, x2, ...xn where xi can be 0 or 1. I want to calculate the following:
sum{ beta[a+sum(xi), b+n-sum(xi) ]* [ (1-x1)dnorm(0,1)+x1dnorm(2,1) ]* [
(1-x2)dnorm(0,1)+x2dnorm(2,1) ]* ...* [ (1-xn)dnorm(0,1)+xndnorm(2,1) ] }
The problem I have is that you're taking the product of the n normal
mixtures and then
you somehow want to sum over them; if
this is meant to be a likelihood function, then it seems that you would
want
beta(a+sum(x), b+n-sum(x)) * prod((1 - x) * dnorm(x, 0, 1) + x * dnorm(x,
2, 1))
instead. As I said up front, it's not clear what you really want, so I
leave you with
multiple possibilities in the hope that one of them is what you need...
The sum in the beginning is over all 2^n possible values for the vector
x1,
x2, ...xn .
Thank you very much!
Hannah
[[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<http://www.r-project.org/posting-guide.html>
and provide commented, minimal, self-contained, reproducible code.
[[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.
Charles C. Berry (858) 534-2098
Dept of Family/Preventive Medicine
E mailto:cbe...@tajo.ucsd.edu UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901
______________________________________________
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.
