Dear Martin, Thank you for the solutions.
I've tried as per your codes and my problem solved :) Thanks a lot! Best regards, Isaudin On Tue, Aug 23, 2016 at 9:06 AM, Martin Maechler <maech...@stat.math.ethz.ch > wrote: > >>>>> Isaudin Ismail <isau...@gmail.com> > >>>>> on Thu, 18 Aug 2016 17:03:50 +0100 writes: > > > Dear Martin, Following my earlier question on "error while > > fitting gumbel copula", I have also crated a new gist at > > https://gist.github.com/anonymous/0bb8aba7adee550d40b840a47d8b7e25 > > for easy checking and copying codes. > > > I got no problem fitting other Archimedean copulas except > > gumbel copula as per my code I used above. > > > Appreciate your kind help. > > > Many thanks, Isaudin > > > > On Mon, Aug 15, 2016 at 4:28 PM, Isaudin Ismail > > <isau...@gmail.com> wrote: > > >> Dear Dr. Martin, > >> > >> I'm glad that you replied to my queries. > >> > >> As advised, I have prepared the following: > > MM: I'm including (cut'n'pasting) my commented and augmented version here: > > ------------------------------------------------------------------ > > ## From: Isaudin Ismail <isau...@gmail.com> > ## To: Martin Maechler <maech...@stat.math.ethz.ch> > ## CC: <r-help@r-project.org> > ## Subject: Re: [R] Error while fitting gumbel copula > ## Date: Mon, 15 Aug 2016 16:28:14 +0100 > > ## Dear Dr. Martin, > > ## I'm glad that you replied to my queries. > ## As advised, I have prepared the following: > > library(copula) > > ## 5 series of data, A, B, C, D and E > A <- c(0.849420849, 0.900652985, 0.97144217, 0.817888428, 0.877901578, > 1.070040669, 0.889742431, 0.87588968, 0.853541938, 0.848664688, > 0.876830319, 0.749582638, 0.818515498, 0.890997174, 0.794766966, > 0.784794851, 0.814858959, 1.074396518, 0.83752495, 0.894341116, > 0.880375293, 0.900816803) > > B <- c(0.479850746, 0.652111668, 1.880607815, 0.579902303, 0.50669344, > 0.747560182, 0.701754386, 0.48969697, 0.346751006, 0.379234973, > 0.862691466, 0.328280188, 0.317312661, 0.534438115, 0.487002653, > 0.335043612, 0.373346897, 0.627520161, 0.792114695, 0.938253012, > 0.444553967, 0.625972763) > > C <- c(0.693491124, 0.866523143, 4.585714286, 1.512055109, 0.387755102, > 0.513435701, 0.76252505, -0.113113113, 0.338521401, 0.333951763, > 0.668755595, 0.401273885, 0.419868791, 0.272885789, 0.541541542, > 0.32751938, 0.386409736, 0.957446809, 0.861195542, 1.531632653, > 0.431610942, 1.226470588) > > D <- c(0.807792208, 0.548547718, 0.738232865, 0.542247744, 1.088964927, > 0.862385321, 0.60720268, 1.000816993, 0.699289661, 0.41723356, > 0.604037267, 0.605003791, 0.698940998, 0.764792899, 0.647897898, > 0.825256975, 0.767476085, 0.941391941, 0.889547813, 0.324503311, > 0.942435424, 0.740686633) > > E <- c(1.077598829, 0.318507891, 1.152616279, 0.930397727, 1.515994437, > 0.940689655, 0.880886427, 1.054274084, 1.067282322, 0.677419355, > 0.966233766, 0.761029412, 1.05734767, 0.615925059, 1.061988304, > 1.07184241, 1.058890147, 1.123873874, 1.304891923, -0.069584736, > 1.172757475, 0.501096491) > > require(copula) > gumbel.copula <- gumbelCopula(dim = 2) > p <- pobs(cbind(D + E, A + B+ C )) > > fit.gumbel <- fitCopula(gumbel.copula, p, method = "ml") > > ## The error is here when trying to fit the gumbel copula > > # I got the following error: > ## Error in optim(start, loglikCopula, lower = lower, upper = upper, > method = > ## method, : > ## non-finite finite-difference value [1] > ## In addition: Warning message: > ## In .local(copula, tau, ...) : tau is out of the range [0, > 1] > > ## MM: my version of copula gives the error message "some tau < 0" > ## -- --------- > ## and indeed: > (tau.p <- cor(p[,1], p[,2], method="kendall")) > ## [1] -0.1428571 > ## ^--------- Kendall's tau is = - 1/7 < 0 ... and that is not > good for Gumbel! > > plot(p) > > ##---------------------------------------------------------- > ------------------- > > So, you tried fitting to *negatively* correlated data, and if > you use the default instead of "ml" the copula is fit, and uses > param = 1 (which corresponds to the *independence* copula: > Because among all the (weakly) positively correlated gumbel > copulas, the boundary case, param = 1 (<==> tau = 0) is the best > fitting. > > What you can do is to "rotate" the data (actually mirror it), > and fit a gumbel copula, which now works nice and easily : > > p2 <- p; p2[,2] <- 1-p2[,2] > (tau.p2 <- cor(p2, method="kendall")) > ## --> now positively correlated > ## ---> > gumb.ml.p2 <- fitCopula(gumbel.copula, p2, method = "ml") > summary(gumb.ml.p2) # looks fine now : > > Call: fitCopula(copula, data = data, method = "ml") > Fit based on "maximum likelihood" and 22 2-dimensional observations. > Gumbel copula, dim. d = 2 > Estimate Std. Error > param 1.121 0.209 > The maximized loglikelihood is 0.1839 > Optimization converged > Number of loglikelihood evaluations: > function gradient > 6 6 > > > --- > > The next version of copula --- or the R-forge one, if you are > interested will support fitting "rotated" copulas, > so you fit a rotated gumbel copula, flipping, the 2nd coordinate > (but not the first): > > rotGcop <- rotCopula(gumbelCopula(), flip=c(FALSE,TRUE)) > > f2 <- fitCopula(rotGcop, data = p) # default method: "mlp" > summary(f2) > > ## Call: fitCopula(copula, data = data) > ## Fit based on "maximum pseudo-likelihood" and 22 2-dimensional > observations. > ## rotCopula copula: Gumbel copula, dim. d = 2 > ## Estimate Std. Error > ## param 1.121 0.225 > ## The maximized loglikelihood is 0.1839 > ## Optimization converged > ## Number of loglikelihood evaluations: > ## function gradient > ## 6 6 > > ------- > > Best regards, > > Martin Maechler, > ETH Zurich > [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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.
