[R] how to draw a multivariate function

[Deepayan Sarkar]

On Saturday 16 October 2004 23:12, Sun wrote:
> Hi, All:
>
> Thanks. Here is the code
>
> n = 30
> lamdaa = 4
> lamdab = 1.5
>
> pa = lamdaa/n
> pb = lamdab/n
>
> x <- seq(0, n/2, len = n/2+1)
> y <- seq(0, n/2, len = n/2+1)

Have you looked at what values of x and y these produce? They include 
non-integer values. Are you sure you want that?

> f = factorial(n)/ (factorial(x) * factorial(y) * factorial (n-x-y))*
> pa^x * pb^y * ((1-pa-pb)^(n-x-y))
> wireframe(f ~ x * y, shade = TRUE)
>
> The above cannot show anything.
> Just le t you know that now I changed to cloud, it can display
> something :) cloud(f ~ x * y, shade = TRUE)
>
> I have questions:
>
> 1.
> what does x*y mean here? I don't think it is a vector dot
> multiplication. I guess it will creat all rows of x and y for all
> possible combinations? Why wireframe cannot show here?

Why guess instead of reading the documentation and looking at the 
examples? There's a very relevant example in the help page for 
wireframe.

You clearly want to evaluate 'f' at all combinations of x and y, yet you 
seem to be evaluating it only along the diagonal (x = y). The correct 
way to do this is (as studying the examples should have suggested to 
you):

g <- expand.grid(x = seq(0, n), y = seq(0, n))
g$z <- dtrinom(g$x, g$y)
wireframe(z ~ x * y, data = g, shade = TRUE)

where dtrinom could be defined as 

dtrinom <- function(x, y)
{
    ifelse(x + y > n,
           NA,
           factorial(n)/ (factorial(x) * 
           factorial(y) * factorial (n-x-y))*
           pa^x * pb^y * ((1-pa-pb)^(n-x-y)))
}

although I would suggest working on the log scale for numerical 
stability:

dtrinom <- function(x, y)
{
    ifelse(x + y > n,
           NA,
           exp((lfactorial(n) - lfactorial(x) -
                lfactorial(y) - lfactorial(n-x-y)) + 
               x * log(pa) +  y * log(pb) +
               (n-x-y) * log(1-pa-pb)))
}


> 2.
> How to show the value on the cloud plot? I have no idea of how much
> the data value is from the plot.

Read the documentation for the 'scales' argument.

> 3. Where can I get resources of R? The help file seems not very
> helpful to me. For example, the lm () function, its weighted least
> square option does not say clearly the weight = standard deviation.
> It said it is to minimize sum w*error^2, which mislead us to think it
> takes variance. I have to ask experienced people. And everytime the
> answer depends on luck.

It's too bad you feel that way. Statistics, and in particular linear 
modeling, is a non-trivial subject, and R documentation is not supposed 
to serve as a textbook. If you don't understand what "minimizing 
'sum(w*e^2)'" means, you really do need help from 'experienced people'. 
Alternatively, look at the references listed in the help page for lm.

Hope that helps,

Deepayan


>
> Thanks,
>
> ----- Original Message -----
> From: "Andrew Ward" <Andrew.Ward at qsa.qld.edu.au>
> To: "Sun" <sun at cae.wisc.edu>
> Sent: Saturday, October 16, 2004 10:15 PM
> Subject: RE: [R] how to draw a multivariate function
>
> > Dear Sun,
> >
> > Could you please provide an example that can be run
> > by readers of the list? What you've given is
> > missing at least n and pa.