[R] Data Gaps
R community, I am trying to write a code that fills in data gaps in a time series. I have no R or statistics background at all but the use of R is proving to be a large portion of my PhD research. So far my code identifies where and the number of new entries required but I do not know how to add additional rows or columns into an array. Any advice on how this can be done? Here is an example: H [0.88 0.72 0.89 0.93 1.23 0.86] T [7.14 7.14 7.49 8.14 7.14 7.32] O [0 0 0 2 0 0] This says that in order to complete the data set 2 entries are required prior to H[4] and T[4] i.e. where O = 2. Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/Data-Gaps-tp2993317p2993317.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] Date Time Objects
I am trying to convert an array from numeric values back to date and time
format. The code I have used is as follows;
for (i in 0:(length(DateTime3)-1)) {
DateTime3[i] <- (strptime(start, "%m/%d/%Y %H:%M")+ i*interval)
where start <- [1] "1/1/1981 00:00"
However the created array (DateTime3) contains [1,] 347156400 347157600
347158800 34716 347161200 347162400 347163600 NA
Does anyone know how I can change DateTime 3 to the same format as start?
Thanks,
Doug
--
View this message in context:
http://r.789695.n4.nabble.com/Date-Time-Objects-tp2993524p2993524.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@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.
Re: [R] Data Gaps
Dennis,
Thanks for that. The problem now is that I am trying to use it in a for
loop. Based on the example before, 2 entries are required after H[3] as
specified by O. The problem is that when inserting values the length of the
array changes so I don't know how to tell the loop that. This is what I
have:
for (i in 1:length(H)) {
if (o[i]>0) {append(H, -1, after = i-1)}
}
Any ideas how to create a loop that allows the inclusion of more than 1
value?
Cheers,
Doug
djmuseR wrote:
>
> Hi:
>
> Perhaps
>
> ?append
>
> for simple insertions...
>
> HTH,
> Dennis
>
> On Wed, Oct 13, 2010 at 1:24 AM, dpender wrote:
>
>>
>> R community,
>>
>> I am trying to write a code that fills in data gaps in a time series. I
>> have no R or statistics background at all but the use of R is proving to
>> be
>> a large portion of my PhD research.
>>
>> So far my code identifies where and the number of new entries required
>> but
>> I
>> do not know how to add additional rows or columns into an array. Any
>> advice
>> on how this can be done?
>>
>> Here is an example:
>>
>> H [0.88 0.72 0.89 0.93 1.23 0.86]
>> T [7.14 7.14 7.49 8.14 7.14 7.32]
>> O [0 0 0 2 0 0]
>>
>> This says that in order to complete the data set 2 entries are required
>> prior to H[4] and T[4] i.e. where O = 2.
>>
>> Thanks,
>>
>> Doug
>> --
>> View this message in context:
>> http://r.789695.n4.nabble.com/Data-Gaps-tp2993317p2993317.html
>> Sent from the R help mailing list archive at Nabble.com.
>>
>> __
>> R-help@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.
>>
>
> [[alternative HTML version deleted]]
>
> __
> R-help@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.
>
>
--
View this message in context:
http://r.789695.n4.nabble.com/Data-Gaps-tp2993317p2993582.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@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.
Re: [R] Date Time Objects
Thanks Henrique, Do you have any idea why the first entry doesn't have the time as the start specified is "1/1/1981 00:00"? Is it something to do with being at midnight? Doug -- View this message in context: http://r.789695.n4.nabble.com/Date-Time-Objects-tp2993524p2993640.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Date Time Objects
Perfect. Thanks -- View this message in context: http://r.789695.n4.nabble.com/Date-Time-Objects-tp2993524p2993663.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Data Gaps
Thanks Dennis. One more thing if you don't mind. How to I abstract the individual H and T “arrays” from f(m,o,l) so as I can combine them with a date/time array and write to a file? Sorry if it’s a simple question but I’m completely new to R. Cheers, Doug -- View this message in context: http://r.789695.n4.nabble.com/Data-Gaps-tp2993317p2994997.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Building ragged dataframe (was Re: Data Gaps)
The H and T values come form the same recording so the gaps that require filling in H are the same as those in T. Therefore the lists will never be of unequal lengths and I can abstract the list to an array using array(unlist(H)) which is a lot simpler. My fault for not clarifying! Finally, the actually dataset that needs filled will have multiple gaps like the one below: H <- c(0.88, 0.72, 0.89, 0.93, 1.23, 0.86, 0.98, 0.85, 1.23) T <- c(7.14, 7.14, 7.49, 8.14, 7.14, 7.32, 7.14, 7.14, 7.14) O <- c(0, 0, 0, 2, 0, 0, 0, 2, 0) Can you advise how to create multiple input lists from this so that I can incorporate that into the function that you defined? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/Data-Gaps-tp2993317p2995263.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] Create Arrays
Hi, For this example: O <- c(0 0 0 2 0 0 2 0) I want to create an array every time O[i] > 0. The array should be in the form; R[j] <- array(-1, dim=c(2,O[i])) i.e. if O[i] > 0 4 times I want 4 R arrays. Does anyone have any suggestions? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/Create-Arrays-tp2996706p2996706.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Create Arrays
Barry, Gerrit,
That was what I am after but unfortunately only the starting point. I am
now trying to amend a function that inserts the R matrices into a dataset in
the correct places:
i.e.
H <- c(0.88, 0.72, 0.89, 0.93, 1.23, 0.86, 0.98, 0.85, 1.23)
T <- c(7.14, 7.14, 7.49, 8.14, 7.14, 7.32, 7.14, 7.14, 7.14)
O <- c(0, 0, 0, 3, 0, 0, 0, 2, 0)
R <- lapply(O[O!=0], function(x){matrix(-1,1,x)})
m <- rbind(H, T)
O2 <- rbind(O,O)
f <- function(x, y, z) {
if(nrow(x) != nrow(y) || nrow(y) != length(z))
stop('unequal numbers of rows among inputs')
out <- matrix(NA,2,length(H)+sum(O))
for(i in 1:2)
out[[i,]] <- append(x[i,], as.numeric(z[[i]]), after = which(y[i,] >
0) - 1)
out
}
f(m, O2, R)
f is the function that requires amendment and was originally written for a
single R value and treated each row of out separately.
I now want to insert R[i] into m before every point that O[i] > 0
Hope this makes sense.
Doug
I don
--
View this message in context:
http://r.789695.n4.nabble.com/Create-Arrays-tp2996706p2996776.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@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.
Re: [R] Create Arrays
Hi Gerrit,
Almost it but I need to insert M[,i] as well as (matrix( -1, nrow( M),
CN[i]) when CN[i] = 0
I know this is not correct but can something like the following be done?
HH <- c(0.88, 0.72, 0.89, 0.93, 1.23, 0.86, 0.98, 0.85, 1.23)
TT <- c(7.14, 7.14, 7.49, 8.14, 7.14, 7.32, 7.14, 7.14, 7.14)
c <- c(0, 0, 0, 2, 0, 0, 0, 2, 0)
TMP <- lapply( seq(c),
function( i, CN, M) {
if( CN[i] == 0) as.matrix( M[, i]) else
(matrix( -1, nrow( M), CN[i]) && as.matrix( M[, i]))
}, CN = c, M = rbind( HH, TT))
do.call( cbind, TMP)
Doug
--
View this message in context:
http://r.789695.n4.nabble.com/Create-Arrays-tp2996706p2997060.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@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.
[R] difftime error
R community, I am trying to create an array of the time differences between datapoints for a very large set. For some reason for 4 of the values the difference has been calculated as NA. Looking at the individual points two of them are "1981-03-29 01:40:00" and "1981-03-29 02:00:00" This is the exact same format as the other points that return values of 20 mins as expected. In such a large dataset (600,000+ entries) I don't know why only 4 would be doing this. Apologies for such a vague question but can anyone shed any light on what may be wrong? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/difftime-error-tp3010105p3010105.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] difftime error
That's it exactly. Do you know how to specify the timezone for Sydney, Australia? -- View this message in context: http://r.789695.n4.nabble.com/difftime-error-tp3010105p3010175.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] datetime objects
I am trying to create an array of date time objects using the strptime function. The first entry is in the form "1/1/1981 0:00" and U use; strptime(datetime, "%m/%d/%Y %H:%M") which gives "1981-01-01 EST". How do I alter this to give me "1981-01-01 0:00 EST"? Thanks in advance, Doug -- View this message in context: http://r.789695.n4.nabble.com/datetime-objects-tp3013417p3013417.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] Clustering
I am looking to use R in order to determine the number of extreme events for a high frequency (20 minutes) dataset of wave heights that spans 25 years (657,432) data points. I require the number, spacing and duration of the extreme events as an output. I have briefly used the clusters function in evd package. Can anyone suggest a more appropriate package to use for such a large dataset? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/Clustering-tp3017056p3017056.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Clustering
That's helpful but the reason I'm using clusters in evd is that I need to specify a time condition to ensure independence. I therefore have an output in the form Cluster[[i]][j-k] where i is the cluster number and j-k is the range of values above the threshold taking account of the time condition. >From this I can get durations easily enough but the spacing is proving quite difficult. The data is for ocean waves and therefore it may be possible that the wave height drops below the threshold for a short period but should still be considered part of the same event, hence the time conditon. Hope this clarifies the problem. Doug -- View this message in context: http://r.789695.n4.nabble.com/Clustering-tp3017056p3018744.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Clustering
Apologies for being vague, The structure of the output is as follows: $ cluster1 : Named num [1:131] 3.05 2.71 3.26 2.91 2.88 3.11 3.21 -1 2.97 3.39 ... ..- attr(*, "names")= chr [1:131] "6667" "6668" "6669" "6670" ... With 613 clusters. What I require is abstracting the first and last value of - attr(*, "names")= chr [1:131] "6667" "6668" "6669" "6670" This will give a start and end point of the cluster and allow for the spacing to be determined. Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/Clustering-tp3017056p3019323.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Clustering
David Winsemius wrote: > > > On Oct 29, 2010, at 12:08 PM, David Winsemius wrote: > >> >> On Oct 29, 2010, at 11:37 AM, dpender wrote: >> >>> Apologies for being vague, >>> >>> The structure of the output is as follows: >> >> Still no code? >> > > I am using the Clusters function from the evd package > >>> >>> $ cluster1 : Named num [1:131] 3.05 2.71 3.26 2.91 2.88 3.11 3.21 >>> -1 2.97 >>> 3.39 ... >>> ..- attr(*, "names")= chr [1:131] "6667" "6668" "6669" "6670" ... >>> >>> With 613 clusters. What I require is abstracting the first and >>> last value >>> of >>> >>> - attr(*, "names")= chr [1:131] "6667" "6668" "6669" "6670" >> >> Those values are in an attribute: > > Corrections: >> >> ? attribute > > ?attributes > >> ? attr >> >> Your specific request may (perhaps) be addressed by something like: >> >> attrnames <- attr(objname["cluster1"], "names") > ^ ^ should be doubled square- > > This works to abstract the part that I am looking for but in order to loop > this over every cluster I need an output object of the same form as > clusters to write the names to. > > brackets >> attrnames[c(1, length(attrnames)] > ^ missing right-paren > > Might work: > attrnames <- attr(clusobj[["cluster1"]], "names") > attrnames[c(1, length(attrnames))] > -- > > David Winsemius, MD > West Hartford, CT > > __ > R-help@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. > > Additionally I can get the output as a matrix in form atomic [1:613] 3.01 4.1 3.04 3.81 3.55 3.37 3.09 4.1 3.61 6.36 ... - attr(*, "acs")= num 47.6 where "acs" is the average size. Each height value in the vector has a corresponding number relating to the location in the dataset. When I change the vector to matrix this looks like the row name but it isn't as rownames(clusters) yields NULL. Do you have any idea how to abstract these values? Doug -- View this message in context: http://r.789695.n4.nabble.com/Clustering-tp3017056p3020216.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] Plot Axes
R community, I am creating a bivariate return level plot by adding calculated return period values as lines onto an existing plot using the following code with the points representing the return periods. plot(H2,D2,pch="+",axes=TRUE) points(H.10,D.10, type="l",col="blue") points(H.20,D.20, type="l",col="green") points(H.50,D.50, type="l",col="red") points(H.100,D.100, type="l",col="orange") The problem is that my return period values are greater than the data values and therefore are partially cut out of the plot. How can I increase the axes limits in order to include all of the return period lines? I've tried ## xis(2,at=seq(35,max(D.100),by=20)) ## but it doesn't work. Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/Plot-Axes-tp3036571p3036571.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Plot Axes
http://r.789695.n4.nabble.com/file/n3036574/plot.jpeg -- View this message in context: http://r.789695.n4.nabble.com/Plot-Axes-tp3036571p3036574.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
Re: [R] Plot Axes
Thanks Jim, The next issue is how i get the RP lines to cut the axes. The lines either cut the x axis at D=35 or the y axis at H=4.5. How do I set the axes so that the lines start and end on the axes? http://r.789695.n4.nabble.com/file/n3037496/plot.jpeg -- View this message in context: http://r.789695.n4.nabble.com/Plot-Axes-tp3036571p3037496.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] clusters in evd package
Hi, I am using the clusters function in the evd package and was wondering if anyone knew how to set the row names as date/time objects. My data has 1 column of date and time and another of wave height (H). My two options are: 1 - When in import the data using read.table make sure that the row.names is set = 1 giving row names to the H data. 2 - Import H from the data file and set the row names using row.names(H) <- x where x is the data/time data. The problem is when I try to use the cluster function on either I get the following error: > Clusters <- clusters(H, 3.0, r=144, plot=TRUE, cmax=TRUE, row.names=TRUE) Error in xy.coords(x, y, xlabel, ylabel, log) : (list) object cannot be coerced to type 'double' In addition: Warning message: In min(diff(xdata)) : no non-missing arguments to min; returning Inf Does anyone know how to set date/time row names in the cluster function? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/clusters-in-evd-package-tp3041023p3041023.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] clusters in evd package
Hi, I am using the clusters function in the evd package and was wondering if anyone knew how to set the row names as date/time objects. My data has 1 column of date and time and another of wave height (H). My two options are: 1 - When in import the data using read.table make sure that the row.names is set = 1 giving row names to the H data. 2 - Import H from the data file and set the row names using row.names(H) <- x where x is the data/time data. The problem is when I try to use the cluster function on either I get the following error: > Clusters <- clusters(H, 3.0, r=144, plot=TRUE, cmax=TRUE, row.names=TRUE) Error in xy.coords(x, y, xlabel, ylabel, log) : (list) object cannot be coerced to type 'double' In addition: Warning message: In min(diff(xdata)) : no non-missing arguments to min; returning Inf Does anyone know how to set date/time row names in the cluster function? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/clusters-in-evd-package-tp3046505p3046505.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] Storm Clustering using clusters in evd
Hi, I am using the clusters function in the evd package in order to determine storm events from a wave time series. So far I have the code working as I want it for wave height on its own but I would now like to include the period as well. The input data is in the form of: H t H t H t so every height measurement has a corresponding period. The storms are defined when the wave height exceeds a certain value and what I am looking to do is to retain the corresponding periods relating to the wave heights in the cluster. This would essentially result in a cluster with 2 variables. Does any one have any ideas? Thanks in advance, Doug -- View this message in context: http://r.789695.n4.nabble.com/Storm-Clustering-using-clusters-in-evd-tp3229951p3229951.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] uniroot
Hi,
I am using the uniroot function in order to carry out a bivariate Monte
Carlo simulation using the logistics model.
I have defined the function as:
BV.FV <- function(x,y,a,A)
(((x^(-a^-1)+y^(-a^-1))^(a-1))*(y^(a-1/a))*(exp(-((1^(-a^-1)+y^(-a^-1))^a)+y^-1)))-A
and the procedure is as follows:
Randomly generate values of A~(0,1), y0 = -(lnA)^-1
Where: A=Pr{Xhttp://r.789695.n4.nabble.com/uniroot-tp3260090p3260090.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@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.
Re: [R] uniroot
Thanks for your advice. There was an error in the equation that is was copying. Doug -- View this message in context: http://r.789695.n4.nabble.com/uniroot-tp3260090p3264288.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] Function within functions and MLE
Hi,
I am trying to determine the MLE of the following function:
http://r.789695.n4.nabble.com/file/n3317341/untitled.bmp
I have defined both parts of the equation as separate functions and looped
over the t and G values to get summations of each part.
The lamda function has 3 unknowns which I am trying to determine using MLE
bub tin order to try and get the overall function working these have been
assigned values at the moment.
So what I am trying to do is write a function l(t|0) that is the sum of two
functions X and Y which are represented by the summations in the above
equation.
My R code is as follows
te <-
c(11,161,337,684,1075,1639,2021,2325,2715,3052,3711,3928,4417,4808,5527,6489,6832,7641,8074)
G <-
c(106,164,314,350,523,317,281,378,327,578,201,427,338,682,952,337,803,424,0)
Poisson.lik <- function(theta0, theta1, theta2, w, t1, t2) {
X <- array(0, dim=c(1,19))
for (i in 1:19) {
X[i] <- function(theta0, theta1, theta2, w, t1, t2)
(((theta0*w*t2)-(theta1*cos(w*t2))+(theta2*sin(w*t2)))/w) -
(((theta0*w*t1)-(theta1*cos(w*t1))+(theta2*sin(w*t1)))/w)
}
XX <- -sum(X)
Y <- array(0, dim=c(1,19))
for (i in 1:19) {
Y[i] <- function(theta0, theta1, theta2, w, t2) theta0 + theta1*sin(w*t2) +
theta2*cos(w*t2)
}
YY <- sum(log(Y))
LL <- XX+YY
return (-LL)
}
optim(c(0,100), Poisson.lik, w=6.28, t1=te[i], t2=te[i]+G[i])
The code isn't working but I am also not sure if this is the best way to go
about it as I am not sure about how to define functions within functions.
What I want out is the MLE for theta, theta1 and theta1.
Can anyone help? Apologies if the code is confusing!
Thanks,
Doug
--
View this message in context:
http://r.789695.n4.nabble.com/Function-within-functions-and-MLE-tp3317341p3317341.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@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.
[R] mle
Hi, I am looking for some help regarding the use of the mle function. I am trying to get mle for 3 parameters (theta0, theta1 and theta2) that have been defined in the the log-likelihood equation as theta0=theta[1], theta1=theta[2] and theta2=theta[3]. My R code for mle is: mle(Poisson.lik, start=list(theta=c(20,1,1), method="Nelder-Mead", fixed=list(w=w, t1=t1, t2=t2)) But I keep getting the following error, Error in eval(expr, envir, enclos) : argument is missing, with no default I am trying to maximise theta starting at (20,1,1) over a fixed range of w, t1 and t2. Can anyone shed some light as to what is going on? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/mle-tp3318857p3318857.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] mle
Hi, I am looking for some help regarding the use of the mle function. I am trying to get mle for 3 parameters (theta0, theta1 and theta2) that have been defined in the the log-likelihood equation as theta0=theta[1], theta1=theta[2] and theta2=theta[3]. My R code for mle is: mle(Poisson.lik, start=list(theta=c(20,1,1), method="Nelder-Mead", fixed=list(w=w, t1=t1, t2=t2)) But I keep getting the following error, Error in eval(expr, envir, enclos) : argument is missing, with no default I am trying to maximise theta starting at (20,1,1) over a fixed range of w, t1 and t2. Can anyone shed some light as to what is going on? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/mle-tp3318856p3318856.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] MLE
Hi, I am having problems carrying out a mle for 3 parameters in a non-homogenous poisson process. I am trying to use the optim function to minimise the -ve log-likelihood. When I use assumed values of my three parameters (20,1,1) the -ve log-likelihood function returns a value of 1309122 but I i then use these values as a starting point in the optim function the parameter estimates and the function value are much lower. Below is a summary of the output: > optim(c(20,1,1), fn=Poisson.lik, gr=NULL, method="Nelder-Mead",w=w, t1=t1, > t2=t2) $par [1] 0.004487104 -2.657468035 12.186003355 $value [1] 289.6901 $counts function gradient 220 NA $convergence [1] 0 $message NULL There were 50 or more warnings (use warnings() to see the first 50) Where the warnings are: 1: In log(((theta0 * w * t2[i]) - (theta1 * cos(w * t2[i])) + ... : NaNs produced I thought I was using optim correctly but obviously not! Does anyone have any suggestions as to what to try? Thanks, Doug -- View this message in context: http://r.789695.n4.nabble.com/MLE-tp3320852p3320852.html Sent from the R help mailing list archive at Nabble.com. __ R-help@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.
[R] If Statement
Hi,
I am having some problems using the if statement correctly. I have used it
many times previously so I dona't know what is different with this case.
Here is my problem:
I have a 1X10 matrix of values as follows:
> H.MC
[,1]
[1,] 4.257669
[2,] 7.023242
[3,] 4.949857
[4,] 5.107000
[5,] 4.257669
[6,] 4.257669
[7,] 4.257669
[8,] 4.257669
[9,] 4.257669
[10,] 4.257669
What I want to do is replace all the values if 4.257669 with a random number
between 3 and 4.5. To do this I have:
H.MC.fin <- matrix(0,10,1)
for (j in 1:10) {
if(H.MC[j] == 4.257669) H.MC.fin[j] <-runif(1,3,4.5) else H.MC.fin[j] <-
H.MC[j]
}
This doesn't seem to do anything and H.MC.fin is the same as H.MC.
Does anyone know what I am doing wrong?
Thanks,
Doug
--
View this message in context:
http://r.789695.n4.nabble.com/If-Statement-tp3341167p3341167.html
Sent from the R help mailing list archive at Nabble.com.
__
R-help@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.
