[R-sig-eco] rank lognormal

[Mario José]

Hi Oksanen,

thank you for help. I had not thought about estimated parameters. I need search how to do this.

Converting fitted values to rank abundance from same models is very difficult. I can't found how to converting fitted logserie, for instance, to rank abundance. In May (1975) text, he explain how to do and in Stevens (2009 - Primer of ecology with R) there are a implementation in R. But it result in values to plot and not a list of values ranked. I try too implement May recommendation, but I cant do this.

The recommendation is use rank abundance to compare different fitted models with data. My intention is to use ranks of models to test my data. I can not use AIC approach because I am using a third model (Tokeshi) that not have implementations for this.

Thank you again.

Regards,

Mario


----- Mensagem original -----
> De: Jari Oksanen <jari.oksanen at oulu.fi>
> Para: Mario José <mariognu-eco at yahoo.com.br>; "r-sig-ecology at r-project.org" <r-sig-ecology at r-project.org>
> Cc: 
> Enviadas: Quarta-feira, 30 de Maio de 2012 10:45
> Assunto: RE: [R-sig-eco] rank lognormal
> 
> Mario,
> 
> There are two important differences in your code and in Wilson's text:
> 
> 1) \Phi^{-1} (or your fi^-1) is standard normal, i.e, in R qnorm(y, 0, 1) 
> instead of your qnorm(y, mulog, sdlog). The mulog and sdlog come back to the 
> play with your multiplication and addition in the last line.
> 
> 2) the term (S - i + 0.5)/S is the same as reverse ppoints function in R. Now 
> you use reverse as the last term of the multiplication, but non-reversed in your 
> fi() function. In fact f(p[i]) * p[i] should do the trick, as these two p[i] 
> terms must be the same -- one cannot be reversed and another non-reversed.
> 
> A working function with minimal changes is:
> 
> rank.lognormal <-
> function(x){
>   S <- length(x);
>   xlog <- log(x);
>   p <- rev(ppoints(xlog));
>   mulog <- mean(xlog);
>   sdlog <- sd(xlog);
>   fi <- function(y){ qnorm(y)  };
>   sapply(1:S, FUN=function(i){ exp(1) ^ (mulog + sdlog * fi(p[i]) * (S-i+0.5)/S) 
> });
> }
> 
> Then some R things.
> 
> 1) do not use ";" to end the line
> 2) exp(1)^(x) == exp(x) -- use the latter form
> 3) you have unnecessary sapply: you can directly work with vectors without 
> accessing individual vector elements.
> 
> Actually, you write the function as a one-liner, and with some elaboration you 
> can streamline it otherwise, but a direct re-implementation of Bastow's 
> formulation and a re-incarnation of your original code is:
> 
> rank.lognormal <-
> function (x) 
> {
>    p <- rev(ppoints(x))
>    logx <- log(x)
>    exp(mean(logx) + sd(logx) * qnorm(p) * p)
> }
> 
> This only returns the fitted values, like your original code. You need some 
> changes to also return the estimated parameters.
> 
> The fitted totals are now closer to observed totals, but not identical.
> 
> Cheers, Jari Oksanen
> ________________________________________
> From: r-sig-ecology-bounces at r-project.org [r-sig-ecology-bounces at r-project.org] 
> on behalf of Mario José [mariognu-eco at yahoo.com.br]
> Sent: 25 May 2012 16:13
> To: r-sig-ecology at r-project.org
> Subject: [R-sig-eco] rank lognormal
> 
> Hi all,
> 
> I trying develop a code to implement a recommendation of Wilson et al (1991). In 
> they text formula:
> 
> "...
> 
> LnA' = fitted mean ln abundance; sigma = fitted standar deviation of ln 
> abundance.
> 
> In ranked-abundance terms:
> 
> lnAi = lnA' + sigma * fi^-1 * (S-i+0.5)/S
> 
> Where
> S = number of species; Ai = abundance of ith out of S species; fi^-1 =
> inverse cumulative distribution function of a norma distribution, i.e.
> the ln abundance at which the area under the normal curve is the value
> indicated.
> ..."
> 
> I developed this code:
> 
> rank.lognormal<-function(x){
>   S <- length(x);
>   xlog <- log(x);
>   p <- ppoints(xlog);
>   mulog <- mean(xlog);
>   sdlog <- sd(xlog);
>   fi <- function(y){ qnorm(y, mulog, sdlog)  };
>   sapply(1:S, FUN=function(i){ exp(1) ^ (mulog + sdlog * fi(p[i]) * (S-i+0.5)/S) 
> });
> }
> 
> bci<-c(25, 24, 22, 21, 18, 17, 15, 14, 14, 13, 13, 12, 12, 11, 11,
>        10, 9, 8, 7, 6, 6, 6, 6, 5, 5, 5, 5, 5, 5, 4, 4, 4, 4, 4, 4,
>        3, 3, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
>        2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
>        1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1);
> 
> sum(rank.lognormal(bci))
> plot(log(rank.lognormal(bci)))
> 
> But the sum of fitted rank don't match with the sum of rank abundance. The 
> curve don't like lognormal too.
> 
> I know rad.lognormal of vegan package, but I'd like follow Wilson 
> recommendation and understand how this procedure works.
> 
> Thank you.
> 
> Best,
> 
> Mario
> 
> Wilson et al 1991 Methods for fitting dominance/diversity curves - Journal of 
> Vegetation Science.
> http://onlinelibrary.wiley.com/doi/10.2307/3235896/abstract
> 
> ...................................................
> 
> Master student in Ecology
> 
> Institute of Biology, Dept. Plant Biology, Ecology Lab.
> State University of Campinas - UNICAMP
> Campinas, São Paulo, Brazil
> 
> 
> _______________________________________________
> R-sig-ecology mailing list
> R-sig-ecology at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
>