[R] misleading example or ...

mauede at alice.it mauede at alice.it
Fri Feb 20 11:18:08 CET 2009


I am always fascinating by good programming techniques. R contains a lot of very good examples I have been learning from.
Since I am using some functions from package "wmtsa", I though I could borrow the elegant example, in the documentation of function "wavBestBasis" to compute the entropy from the Wavelet Transform coefficients.
In the following I have pasted the excerpt, from "wmtsa" on-line documentation, that implements the entropy calculation:
## define an entropy cost functional 
"entropy" <- function(x){
    iz <- which(x==0)
    z <- -x^2 * log(x^2)
    if (length(iz))
       z[iz] <- 0
    sum(z)
}

To my best recollection, Shannon's entropy operates on probabilities therefore requesting a normalization step.
I have written my simple code that implements Shannon's entropy:
# Shannon's entropy 
# input: vector x
 sum(x^2)
 pi <- (x^2)/sum(x^2)
 pi
 -pi*log(pi)
 -sum(pi*log(pi))

I have  tested both entropy realizations on a simple case:
 >  x<- 1:10
>  x
 [1]  1  2  3  4  5  6  7  8  9 10

# ENTROPY FROM WMTSA EXAMPLE
> entropy(x)
[1] -1552.495

# SHANNON ENTROPY
>  sum(x^2)
[1] 385
>  pi <- (x^2)/sum(x^2)
>  pi
 [1] 0.002597403 0.010389610 0.023376623 0.041558442 0.064935065 0.093506494 0.127272727
 [8] 0.166233766 0.210389610 0.259740260
>  -pi*log(pi)
 [1] 0.01546297 0.04744882 0.08780304 0.13218305 0.17755633 0.22158462 0.26236293 0.29828326
 [9] 0.32795410 0.35014887
>  -sum(pi*log(pi))
[1] 1.920788

Needless to say, when I calculate the entropy on wavelet coefficients to the purpose of selecting the best wavelet  family on the basis of the "entropy" cost function, I get different results depending on which entropy implementation I use ...

Your comments and thoughts are very welcome.














e tutti i telefonini TIM!
Vai su 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: wmtsa.pdf
Type: application/pdf
Size: 551276 bytes
Desc: wmtsa.pdf
URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20090220/d545ae8b/attachment-0002.pdf>


More information about the R-help mailing list