From: Ravi Varadhan <rvaradhan_at_jhmi.edu>

Date: Sat 27 May 2006 - 02:18:07 EST

Ravi Varadhan, Ph.D.

Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University

Ph: (410) 502-2619

Fax: (410) 614-9625

Email: rvaradhan@jhmi.edu

Webpage:http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html

> -----Original Message-----

*> From: r-help-bounces@stat.math.ethz.ch [mailto:r-help-
*

*> bounces@stat.math.ethz.ch] On Behalf Of Chaouch, Aziz
*

*> Sent: Friday, May 26, 2006 12:02 PM
*

*> To: Dimitrios Rizopoulos
*

*> Cc: r-help@stat.math.ethz.ch
*

*> Subject: Re: [R] Maximum likelihood estimate of bivariatevonmises-
*

*> weibulldistribution
*

*>
*

*> Hi,
*

*>
*

*> I'm still strugling with this copula model but this seems to be the way
*

*> to go. I'm now trying to model the marginal distributions and and for
*

*> wind direction I use a mixture of 2 von mises. I'd like to estimate all
*

*> the parameters (m1,m1,kappa1,kappa2,p) by maximizing the likelihood but
*

*> I don't know how to define the likelihood (or log-likelihood) of a
*

*> mixture of 2 Von Mises to use it with the function fitDistr in the MASS
*

*> package. Can you help me define this likelihood function and use it
*

*> through the fitDistr function?
*

*>
*

*> Thanks,
*

*>
*

*> Aziz
*

*>
*

*> -----Original Message-----
*

*> From: Dimitrios Rizopoulos [mailto:Dimitris.Rizopoulos@med.kuleuven.be]
*

*> Sent: May 12, 2006 4:35 PM
*

*> To: Chaouch, Aziz
*

*> Subject: RE: [R] Maximum likelihood estimate of bivariate
*

*> vonmises-weibulldistribution
*

*>
*

*> look at the following code:
*

*>
*

*> library(copula)
*

*> par(mfrow = c(2, 2))
*

*> x <- mvdc(normalCopula(sin(0.5 * pi /2)), c("norm", "norm"),
*

*> list(list(mean = 0, sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc,
*

*> xlim = c(-2.7, 2.7), ylim = c(-2.7, 2.7))
*

*>
*

*> x <- mvdc(frankCopula(5.736276), c("norm", "norm"), list(list(mean = 0,
*

*> sd = 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7),
*

*> ylim = c(-2.7, 2.7))
*

*>
*

*> x <- mvdc(gumbelCopula(2), c("norm", "norm"), list(list(mean = 0, sd =
*

*> 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim
*

*> = c(-2.7, 2.7))
*

*>
*

*> x <- mvdc(claytonCopula(2), c("norm", "norm"), list(list(mean = 0, sd =
*

*> 1), list(mean = 0, sd = 1))) contour(x, dmvdc, xlim = c(-2.7, 2.7), ylim
*

*> = c(-2.7, 2.7))
*

*>
*

*>
*

*> the values of the association parameter I've chosen in each copula
*

*> correspond to Kendall's tau 0.5; assuming also standard normal
*

*> marginal distributions look at the different shapes you get!
*

*>
*

*> If possible try something similar for you case (i.e., using von Mises
*

*> and Weibull marginals) and check if the association shape for a
*

*> specific copula is more appropriate for your application. If this is
*

*> not possible fit models assumig different copulas and check which one
*

*> provides a better fit to your data.
*

*>
*

*> I hope it helps.
*

*>
*

*> Best,
*

*> Dimitris
*

*>
*

*>
*

*>
*

*> Quoting "Chaouch, Aziz" <achaouch@NRCan.gc.ca>:
*

*>
*

*> > Hi Dimitris,
*

*> >
*

*> > I'm not sure to understand your suggestion. How would you build that
*

*> > contour plot for a particular copula and on what is computed the
*

*> > Kendall's tau?
*

*> >
*

*> > Thanks,
*

*> >
*

*> > Aziz
*

*> >
*

*> > -----Original Message-----
*

*> > From: Dimitris Rizopoulos
*

*> > [mailto:dimitris.rizopoulos@med.kuleuven.be]
*

*> > Sent: May 12, 2006 9:57 AM
*

*> > To: Chaouch, Aziz; hydinghua@gmail.com
*

*> > Cc: r-help@stat.math.ethz.ch
*

*> > Subject: Re: [R] Maximum likelihood estimate of bivariate
*

*> > vonmises-weibulldistribution
*

*> >
*

*> > the choice of the copula is, in fact, a model selection problem.
*

*> > First, you could have a look at the contour plots of different
*

*> > copulas
*

*> > (preferably for the same value of Kendall's tau), and decide if some
*

*> > of
*

*> > them assume a more appropriate association structure for your
*

*> > application, compared to the others. Afterwards, you may fit various
*

*> > copula functions, check the fit on the data, compute AIC (since
*

*> > these
*

*> > are typically not nested models), etc.
*

*> >
*

*> > regarding the Von Mises distribution and if could be used in mvdc(),
*

*> > that I don't know. It'd be better to contact the copula package
*

*> > maintainer and ask.
*

*> >
*

*> > I hope it helps.
*

*> >
*

*> > Best,
*

*> > Dimitirs
*

*> >
*

*> > ----
*

*> > Dimitris Rizopoulos
*

*> > Ph.D. Student
*

*> > Biostatistical Centre
*

*> > School of Public Health
*

*> > Catholic University of Leuven
*

*> >
*

*> > Address: Kapucijnenvoer 35, Leuven, Belgium
*

*> > Tel: +32/(0)16/336899
*

*> > Fax: +32/(0)16/337015
*

*> > Web: http://www.med.kuleuven.be/biostat/
*

*> > http://www.student.kuleuven.be/~m0390867/dimitris.htm
*

*> >
*

*> >
*

*> > ----- Original Message -----
*

*> > From: "Chaouch, Aziz" <achaouch@NRCan.gc.ca>
*

*> > To: "Dimitris Rizopoulos" <dimitris.rizopoulos@med.kuleuven.be>;
*

*> > <hydinghua@gmail.com>
*

*> > Cc: <r-help@stat.math.ethz.ch>
*

*> > Sent: Friday, May 12, 2006 3:13 PM
*

*> > Subject: RE: [R] Maximum likelihood estimate of bivariate
*

*> > vonmises-weibulldistribution
*

*> >
*

*> >
*

*> > Thanks a lot! I wasn't aware of that copula package and it could well
*

*> > be
*

*> > appropriate to use it for my application. However if I read the
*

*> > copula
*

*> > help correctly, I still need to know what kind of copula to use to
*

*> > link
*

*> > the distribution of wind speeds and directions. Is there a way to
*

*> > determine this in R?
*

*> >
*

*> > Moreover can I use the Von Mises distribution from the circular or
*

*> > CircStats package into the mvdc function of the copula package or
*

*> > does
*

*> > the mvdc function only recognize distributions available "natively"
*

*> > within R?
*

*> >
*

*> > Thanks again to all, your help is highly appreciated for a newbie
*

*> > like
*

*> > me!
*

*> >
*

*> > Regards,
*

*> >
*

*> > Aziz
*

*> >
*

*> > -----Original Message-----
*

*> > From: Dimitris Rizopoulos
*

*> > [mailto:dimitris.rizopoulos@med.kuleuven.be]
*

*> > Sent: May 12, 2006 3:01 AM
*

*> > To: Philip He; Chaouch, Aziz
*

*> > Cc: r-help@stat.math.ethz.ch
*

*> > Subject: Re: [R] Maximum likelihood estimate of bivariate
*

*> > vonmises-weibulldistribution
*

*> >
*

*> >
*

*> > ----- Original Message -----
*

*> > From: "Philip He" <hydinghua@gmail.com>
*

*> > To: "Chaouch, Aziz" <achaouch@nrcan.gc.ca>
*

*> > Cc: <r-help@stat.math.ethz.ch>
*

*> > Sent: Thursday, May 11, 2006 11:21 PM
*

*> > Subject: Re: [R] Maximum likelihood estimate of bivariate
*

*> > vonmises-weibulldistribution
*

*> >
*

*> >
*

*> > > On 5/11/06, Chaouch, Aziz <achaouch@nrcan.gc.ca> wrote:
*

*> > >>
*

*> > >> Hi,
*

*> > >>
*

*> > >> I'm dealing with wind data and I'd like to model their
*

*> > distribution
*

*> > >> in order to simulate data to fill-in missing values. Wind
*

*> > direction
*

*> > >> are typically following a vonmises distribution and wind speeds
*

*> > >> follow a weibull distribution. I'd like to build a joint
*

*> > distribution
*

*> >
*

*> > >> of directions and speeds as a VonMises-Weibull bivariate
*

*> > >> distribution.
*

*> > >
*

*> > >
*

*> > > In order to built a bivariate distribution from two marginal
*

*> > > distributions (wind direction, wind speed) , more information is
*

*> > > needed to specify the relation between these two marginal
*

*> > > distributions.For example, a conditional distribution may help.
*

*> > >
*

*> >
*

*> >
*

*> > An alternative in such cases (i.e., when marginals are available but
*

*> > the
*

*> > joint is difficult to postulate) is to use copulas, which can
*

*> > construct
*

*> > multivariate distributions from univariate marginals. If this is
*

*> > appropriate for this application, the "copula" package might be of
*

*> > help.
*

*> >
*

*> > Best,
*

*> > Dimitris
*

*> >
*

*> > ---
*

*> > Dimitris Rizopoulos
*

*> > Ph.D. Student
*

*> > Biostatistical Centre
*

*> > School of Public Health
*

*> > Catholic University of Leuven
*

*> >
*

*> > Address: Kapucijnenvoer 35, Leuven, Belgium
*

*> > Tel: +32/(0)16/336899
*

*> > Fax: +32/(0)16/337015
*

*> > Web: http://www.med.kuleuven.be/biostat/
*

*> > http://www.student.kuleuven.be/~m0390867/dimitris.htm
*

*> >
*

*> >
*

*> > Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
*

*> >
*

*> >
*

*> > Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
*

*> >
*

*> >
*

*>
*

*>
*

*> Disclaimer: http://www.kuleuven.be/cwis/email_disclaimer.htm
*

*>
*

*> ______________________________________________
*

*> R-help@stat.math.ethz.ch mailing list
*

*> https://stat.ethz.ch/mailman/listinfo/r-help
*

*> PLEASE do read the posting guide! http://www.R-project.org/posting-
*

*> guide.html
*

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sat May 27 02:23:24 2006

Date: Sat 27 May 2006 - 02:18:07 EST

Hi Aziz,

I am attaching a file that contains the functions that I wrote to compute the MLE for a binary vonMises mixture, using the EM algorithm. It also contains a simulation example. Let me know if you have any trouble.

Hope this is helpful,

Ravi.

Ravi Varadhan, Ph.D.

Assistant Professor, The Center on Aging and Health Division of Geriatric Medicine and Gerontology Johns Hopkins University

Ph: (410) 502-2619

Fax: (410) 614-9625

Email: rvaradhan@jhmi.edu

Webpage:http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html

> -----Original Message-----

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sat May 27 02:23:24 2006

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.1.8, at Sat 27 May 2006 - 06:10:28 EST.

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