[R] about the integrate

[Rolf Turner]

Think about what you are actually getting from pnorm().

For each scalar "s" you want to get

     pnorm(s,mu[2],sigma[2])    pnorm(s,mu[3],sigma[3]) 
pnorm(s,mu[4],sigma[4])

(and then take products of things).

But when "s" is a vector you get

pnorm(s[1],mu[2],sigma[2])    pnorm(s[2],mu[3],sigma[3]) 
pnorm(s[3],mu[4],sigma[4])

pnorm(s[4],mu[2],sigma[2])    pnorm(s[5],mu[3],sigma[3]) 
pnorm(s[6],mu[4],sigma[4]) ...

and the product will be nothing like what you want.

Vectors get re-cycled in R.

Note that you get a single scalar quantity from the "prod(1 - 
pnorm(...))" component of
test2().  That single scalar multiplies each entry of the vector 
dnorm(s,mu[1],sigma[1]).

Whereas you want dnorm(s[i],mu[1],sigma[1]) to be multiplied by a 
product involving
only pnorm() terms evaluated at s[i].  That is what Vectorize() arranges 
for you.

If you are going to use R you should learn its syntax and know about 
things like the
re-cycling of vectors.

     cheers,

     Rolf Turner

On 11/22/13 04:46, dan wang wrote:
> Hi all,
>
> I tried below two methods to calculate the integrate of a same function.
> Two different results are given.
> The first one should be the right answer.
> Can any one help to explain why?
> Another issue is, the first one is not convenient as I have to update the
> mu and sigma outside the function. How can I modify the second one so I
> could issue the parameters in the integrate function.
> Thanks in advance,
>
> Dan
>
> test <- function(s){
>         prod(1-pnorm(s,mu[-1],sigma[-1]))*dnorm(s,mu[1],sigma[1])
>
> }
> testV <- Vectorize(test)
> mu=c(0,0,0,0)
> sigma=c(1,1,1,1)
> integrate(testV,lower=-Inf,upper=Inf)$value
>
>
> test2 <- function(s, mu, sigma){
>         prod(1-pnorm(s,mu[-1],sigma[-1]))*dnorm(s,mu[1],sigma[1])
>
> }
> integrate(test2,mu=c(0,0,0,0),sigma=c(1,1,1,1),lower=-Inf,upper=Inf)$value