Re: [R] S-PLUS code in R

From: Fotis Papailias <fpapailias_at_gmail.com>
Date: Sat, 26 Jul 2008 13:28:06 +0100

I made some experiments and I conclude that if I change the interval (e.g. say c(-100,100)!!!) then the minimum value changes in respect to the simulated data.

So, probably the problem lies in the "optimal" values of the interval?

On Sat, Jul 26, 2008 at 1:00 PM, Fotis Papailias <fpapailias_at_gmail.com>wrote:

> Hi,
>
> thanks for your message.
>
> You mean to rewrite the function like that:
>
> lw <- function(d, x, im)
> {
> peri1 <- per(x)
> len <- length(x)
> m <- len/im
> peri <- peri1[2:(m+1)]
> z <- c(1:m)
> freq <- ((2*pi)/len) * z
> result <- log(sum(freq^(2*d-1)*peri))-(2*d)/m * sum(log(freq))
> }
>
>
> and then to use:
>
> target <- function(d) lw(x, d, im=2)
> k <- optimize(target, interval=c(0, 0.5))
>
> or just,
>
> k <- optimize(lw, x, im)
>
> But I still get the same values
>
> $minimum
> [1] 0.4999542
>
> $objective
> [1] 2.739509
>
> after the optimization, no matter my "x" series. It seems that I am still
> doing something wrong in the optimize() or in the result() function inside
> the loop. But comparing it to the original S-PLUS code, it seems like it is
> fine.
>
> Send me if you have any other ideas please.
>
> Thanks again.
>
> fotis
>
>
> On Sat, Jul 26, 2008 at 12:48 PM, Duncan Murdoch <murdoch_at_stats.uwo.ca>wrote:
>
>> On 26/07/2008 7:40 AM, Fotis Papailias wrote:
>>
>>> Dear R-users,
>>>
>>> I have sent another mail some hour ago about a matlab Code I was trying
>>> to
>>> translate in R.
>>>
>>> Actually I have found a simpler code originally written in S-PLUS for the
>>> same function.
>>> Author's page -> http://math.bu.edu/people/murad/methods/locwhitt/
>>>
>>> =============================================================
>>>
>>> rfunc_function(h, len, im, peri)
>>> # h -- Starting H value for minimization.
>>> # len -- Length of time series.
>>> # im -- Use only len/im frequencies.
>>> # peri -- Periodogram of data.
>>> {
>>> m <- len %/% im
>>> peri <- peri[2:(m + 1)]
>>> z <- c(1:m)
>>> freq <- (2 * pi)/len * z
>>> result <- log(sum(freq^(2 * h - 1) * peri)) - (2 * h)/m *
>>> sum(log(freq)
>>> ) # cat("H = ", h, "R = ", result, "\n")
>>> drop(result)
>>> }
>>>
>>>
>>> locwhitt_function(data, h = 0.5, im = 2)
>>> # data -- Time series.
>>> # h -- Starting H value for minimization.
>>> # im -- Use only N/im frequencies where N is length of series.
>>>
>>> {
>>> peri <- per(data)
>>> len <- length(data)
>>> return(nlminb(start = h, obj = rfunc, len = len, im = im, peri =
>>> peri)$
>>> parameters)
>>> }
>>> ===============================================================
>>>
>>> The author who has written the above S-PLUS code uses two functions (with
>>> the locwhitt_function he lets the user to define the data and the
>>> parameters
>>> and with the rfunc_function he does the minimization.)
>>>
>>> Mine translation is in R is:
>>>
>>> where I use a joint function compared to the the above author
>>>
>>>
>>> ================================================================
>>>
>>> lw <- function(x, d, im)
>>> {
>>> peri1 <- per(x)
>>> len <- length(x)
>>> m <- len/im
>>> peri <- peri1[2:(m+1)]
>>> z <- c(1:m)
>>> freq <- ((2*pi)/len) * z
>>> result <- log(sum(freq^(2*d-1)*peri))-(2*d)/m * sum(log(freq))
>>> }
>>>
>>> =================================================================
>>>
>>> which seems to run ok.
>>>
>>> But when I do
>>>
>>> k <- optimize(lw, x, im=2, interval=c(0, 0.5))
>>>
>>> I always get the same result no matter the (simulated) data in x!
>>>
>>> The parameter of interest to be minimized is "d". Does anyone know how to
>>> edit the function "optimize" so it can work properly?
>>>
>>
>> optimize() is fine, but the way you're calling it is not. It optimizes a
>> function over the first argument. So you could rewrite lw to put d first,
>> or write a new function which calls it, e.g.
>>
>> target <- function(d) lw(x, d, im)
>>
>> and then
>>
>> optimize(target, interval=c(0, 0.5))
>>
>> Because target is defined in the global environment, it will look there
>> for x and im, and you don't need to pass them as arguments: unless x and im
>> aren't defined there too!
>>
>> Duncan Murdoch
>>
>
>
>
> --
> fp
>

-- 
fp

	[[alternative HTML version deleted]]

______________________________________________
R-help_at_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.
Received on Sat 26 Jul 2008 - 12:31:51 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Sat 26 Jul 2008 - 13:32:22 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive