[R-sig-phylo] Mismatch distribution (pegas): expected versus empirical and other questions

2013-12-10 Thread Alastair Potts
Hi all,

The mismatch distribution function in pegas (MMD) produces the familiar
histogram and then also a line that I always took as the theoretical
expected distribution under a scenario of constant population size.

I had a quick look at the MMD function and noticed that this line is not
calculated using the theoretical expectations based on theta, but rather is
some sort of moving-average based on the density() function.

Does anyone know how to implement the correct theoretical expectation? I
have dipped into the literature and am now completely lost.

Also, there are a range of statistics (such as the raggedness index) and
Arlequin has somehow managed to calculate a significance value for these
statistics. Has anyone implemented (or can easily implement) these
statistics and significance tests in R?

Thanks in advance for any comments, thoughts and solutions.

Cheers,
Alastair Potts

[[alternative HTML version deleted]]

___
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/


[R-sig-phylo] fitting cuadratic model using gnls

2013-12-10 Thread Jorge Sánchez Gutiérrez
Hello everyone,

 

I’m trying to fit nonlinear models using the function “gnls” in ape. For
example, I run a simple quadratic regression of moult duration (duration)
and latitude (latitude_used) as follows:

 

mod - duration ~ latitude_used +I(latitude_used^2)

init-list(alpha=1,beta=1)

m3 - gnls(mod, data=data, start=init, correlation = dataCor)

 

However, I get the following error message:

 

“Error in gnls(mod, data=data, start = init, correlation = dataCor) : 

  step halving factor reduced below minimum in NLS step”

 

I’ve also tried with several correlation structures and different starting
values for estimation, but without success. I would be very grateful if you
can provide any suggestions.

 

Thanks in advance,

Jorge S. Gutiérrez


[[alternative HTML version deleted]]

___
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/

Re: [R-sig-phylo] fitting cuadratic model using gnls

2013-12-10 Thread Theodore Garland Jr
I don't know, but one thing to watch out for is a high correlation between X 
and Xsquared.
That gives multicollinearity and can lead to numerical problems.
One thing that is commonly done is:
  1.  standardize the X variable (subtract mean and divide by standard 
deviation)
  2.  square that value and use it for the quadratic
This is referred to as an orthogonol polynomial.

Cheers,
Ted

Theodore Garland, Jr., Professor
Department of Biology
University of California, Riverside
Riverside, CA 92521
Office Phone:  (951) 827-3524
Wet Lab Phone:  (951) 827-5724
Dry Lab Phone:  (951) 827-4026
Home Phone:  (951) 328-0820
Skype:  theodoregarland
Facsimile:  (951) 827-4286 = Dept. office (not confidential)
Email:  tgarl...@ucr.edu
http://www.biology.ucr.edu/people/faculty/Garland.html
http://scholar.google.com/citations?hl=enuser=iSSbrhwJ

Inquiry-based Middle School Lesson Plan:
Born to Run: Artificial Selection Lab
http://www.indiana.edu/~ensiweb/lessons/BornToRun.html

Fail Lab: Episode One
http://testtube.com/faillab/zoochosis-episode-one-evolution
http://www.youtube.com/watch?v=c0msBWyTzU0


From: r-sig-phylo-boun...@r-project.org [r-sig-phylo-boun...@r-project.org] on 
behalf of Jorge Sánchez Gutiérrez [jorgesgutier...@unex.es]
Sent: Tuesday, December 10, 2013 9:48 AM
To: r-sig-phylo@r-project.org
Subject: [R-sig-phylo] fitting cuadratic model using gnls

Hello everyone,



I’m trying to fit nonlinear models using the function “gnls” in ape. For
example, I run a simple quadratic regression of moult duration (duration)
and latitude (latitude_used) as follows:



mod - duration ~ latitude_used +I(latitude_used^2)

init-list(alpha=1,beta=1)

m3 - gnls(mod, data=data, start=init, correlation = dataCor)



However, I get the following error message:



“Error in gnls(mod, data=data, start = init, correlation = dataCor) :

  step halving factor reduced below minimum in NLS step”



I’ve also tried with several correlation structures and different starting
values for estimation, but without success. I would be very grateful if you
can provide any suggestions.



Thanks in advance,

Jorge S. Gutiérrez


[[alternative HTML version deleted]]


___
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/


Re: [R-sig-phylo] fitting cuadratic model using gnls

2013-12-10 Thread Ben Bolker
-BEGIN PGP SIGNED MESSAGE-
Hash: SHA1

On 13-12-10 12:57 PM, Theodore Garland Jr wrote:
 I don't know, but one thing to watch out for is a high correlation 
 between X and Xsquared. That gives multicollinearity and can lead 
 to numerical problems. One thing that is commonly done is: 1. 
 standardize the X variable (subtract mean and divide by standard 
 deviation) 2.  square that value and use it for the quadratic This 
 is referred to as an orthogonol polynomial.
 
 Cheers, Ted
 

  Yes.  Therefore, you would be (much) better off wit

 m3 - gls(duration ~ poly(latitude_used,2), data=data,
correlation = dataCor)

(A quadratic model is still a *linear* fit in the statistical sense.)

poly() will generate orthogonal polynomials by default, which will be
more stable and are more appropriate if you're planning to test the
significance of terms.  If it is important to you to have the
parameters interpretable as slope and curvature at the equator (i.e.
where latitude=0), you may be able to get away with

m3 - gls(duration ~ poly(latitude_used,2,raw=TRUE), data=data,
correlation = dataCor)



 Theodore Garland, Jr., Professor Department of Biology University 
 of California, Riverside Riverside, CA 92521 Office Phone:  (951) 
 827-3524 Wet Lab Phone:  (951) 827-5724 Dry Lab Phone:  (951) 
 827-4026 Home Phone:  (951) 328-0820 Skype:  theodoregarland 
 Facsimile:  (951) 827-4286 = Dept. office (not confidential)
 Email: tgarl...@ucr.edu 
 http://www.biology.ucr.edu/people/faculty/Garland.html 
 http://scholar.google.com/citations?hl=enuser=iSSbrhwJ
 
 Inquiry-based Middle School Lesson Plan: Born to Run: Artificial 
 Selection Lab 
 http://www.indiana.edu/~ensiweb/lessons/BornToRun.html
 
 Fail Lab: Episode One 
 http://testtube.com/faillab/zoochosis-episode-one-evolution 
 http://www.youtube.com/watch?v=c0msBWyTzU0
 
  From: 
 r-sig-phylo-boun...@r-project.org 
 [r-sig-phylo-boun...@r-project.org] on behalf of Jorge Sánchez 
 Gutiérrez [jorgesgutier...@unex.es] Sent: Tuesday, December 10, 
 2013 9:48 AM To: r-sig-phylo@r-project.org Subject: [R-sig-phylo] 
 fitting cuadratic model using gnls
 
 Hello everyone,
 
 
 
 I’m trying to fit nonlinear models using the function “gnls” in 
 ape. For example, I run a simple quadratic regression of moult 
 duration (duration) and latitude (latitude_used) as follows:
 
 
 
 mod - duration ~ latitude_used +I(latitude_used^2)
 
 init-list(alpha=1,beta=1)
 
 m3 - gnls(mod, data=data, start=init, correlation = dataCor)
 
 
 
 However, I get the following error message:
 
 
 
 “Error in gnls(mod, data=data, start = init, correlation =
 dataCor) :
 
 step halving factor reduced below minimum in NLS step”
 
 
 
 I’ve also tried with several correlation structures and different 
 starting values for estimation, but without success. I would be 
 very grateful if you can provide any suggestions.
 
 
 
 Thanks in advance,
 
 Jorge S. Gutiérrez
 
 
 [[alternative HTML version deleted]]
 
 
 ___ R-sig-phylo
 mailing list - R-sig-phylo@r-project.org 
 https://stat.ethz.ch/mailman/listinfo/r-sig-phylo Searchable 
 archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/
 

-BEGIN PGP SIGNATURE-
Version: GnuPG v1.4.11 (GNU/Linux)
Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/

iQEcBAEBAgAGBQJSp2DoAAoJEOCV5YRblxUHn4gH/jYtxFo6aj2DlHAPw2FXufA4
CwtJmIHpWlnD5mhg5YqBgJNldn2XPxqx/qEf0mMzRgd2gAyIQiXd91WI6sqkaVhF
FmYSPzClqx14jZSlF64KE54l/8W4mp81rDyYnh+/SSebltCrX8Cv7NqYqvNDhYQC
M7i/avh0NNv/PZuDuUwbj95doQmdeUyeSlig0sxHaP14C7FxzvxODUBhcTB6yQzN
0wuK+iIaUKqSGm3jmcpY08oy0RJ34QfyH175UezRmSaE/RI+XCohhBbdtU5+WgRX
lBRYmwtCIGAJ+DNmhlnhE3u6ZFdnzhw95gYd1LGL3Sr+JZl1U+VIKfIjhlrTgLQ=
=lLum
-END PGP SIGNATURE-

___
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/


Re: [R-sig-phylo] Mismatch distribution (pegas): expected versus empirical and other questions

2013-12-10 Thread David Winter
Hi Alistair,

I  wrote a function to make the standard expectation under constant
population size quite a while ago. It's up here

https://github.com/dwinter/Pegas/blob/new_functions/R/MMD.stable.R

It's been a little while since I've had to think about mismatch
distributions, so, though I'm sure I would checked everything out at
the time, you should probably confirm this function behaves as
expected.

(BTW, Emmanuel, if you want to include any of that function in PEGAS
you are most welcome to)

David Winter




On Tue, Dec 10, 2013 at 6:27 AM, Alastair Potts pott...@gmail.com wrote:
 Hi all,

 The mismatch distribution function in pegas (MMD) produces the familiar
 histogram and then also a line that I always took as the theoretical
 expected distribution under a scenario of constant population size.

 I had a quick look at the MMD function and noticed that this line is not
 calculated using the theoretical expectations based on theta, but rather is
 some sort of moving-average based on the density() function.

 Does anyone know how to implement the correct theoretical expectation? I
 have dipped into the literature and am now completely lost.

 Also, there are a range of statistics (such as the raggedness index) and
 Arlequin has somehow managed to calculate a significance value for these
 statistics. Has anyone implemented (or can easily implement) these
 statistics and significance tests in R?

 Thanks in advance for any comments, thoughts and solutions.

 Cheers,
 Alastair Potts

 [[alternative HTML version deleted]]

 ___
 R-sig-phylo mailing list - R-sig-phylo@r-project.org
 https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
 Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/

___
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/


Re: [R-sig-phylo] fitting cuadratic model using gnls

2013-12-10 Thread Jorge Sánchez Gutiérrez
Thanks Ted and Ben for your useful comments. I used poly() to generate
orthogonal polynomials of degree 2 and it seems to work well. 

Cheers, Jorge

-Mensaje original-
De: Ben Bolker [mailto:bbol...@gmail.com] 
Enviado el: martes, 10 de diciembre de 2013 19:44
Para: r-sig-phylo@r-project.org; Jorge Sánchez Gutiérrez
Asunto: Re: [R-sig-phylo] fitting cuadratic model using gnls

-BEGIN PGP SIGNED MESSAGE-
Hash: SHA1

On 13-12-10 12:57 PM, Theodore Garland Jr wrote:
 I don't know, but one thing to watch out for is a high correlation 
 between X and Xsquared. That gives multicollinearity and can lead to 
 numerical problems. One thing that is commonly done is: 1.
 standardize the X variable (subtract mean and divide by standard
 deviation) 2.  square that value and use it for the quadratic This is 
 referred to as an orthogonol polynomial.
 
 Cheers, Ted
 

  Yes.  Therefore, you would be (much) better off wit

 m3 - gls(duration ~ poly(latitude_used,2), data=data,
correlation = dataCor)

(A quadratic model is still a *linear* fit in the statistical sense.)

poly() will generate orthogonal polynomials by default, which will be more
stable and are more appropriate if you're planning to test the significance
of terms.  If it is important to you to have the parameters interpretable as
slope and curvature at the equator (i.e.
where latitude=0), you may be able to get away with

m3 - gls(duration ~ poly(latitude_used,2,raw=TRUE), data=data,
correlation = dataCor)



 Theodore Garland, Jr., Professor Department of Biology University of 
 California, Riverside Riverside, CA 92521 Office Phone:  (951)
 827-3524 Wet Lab Phone:  (951) 827-5724 Dry Lab Phone:  (951)
 827-4026 Home Phone:  (951) 328-0820 Skype:  theodoregarland
 Facsimile:  (951) 827-4286 = Dept. office (not confidential)
 Email: tgarl...@ucr.edu
 http://www.biology.ucr.edu/people/faculty/Garland.html
 http://scholar.google.com/citations?hl=enuser=iSSbrhwJ
 
 Inquiry-based Middle School Lesson Plan: Born to Run: Artificial 
 Selection Lab
 http://www.indiana.edu/~ensiweb/lessons/BornToRun.html
 
 Fail Lab: Episode One
 http://testtube.com/faillab/zoochosis-episode-one-evolution
 http://www.youtube.com/watch?v=c0msBWyTzU0
 
  From: 
 r-sig-phylo-boun...@r-project.org
 [r-sig-phylo-boun...@r-project.org] on behalf of Jorge Sánchez 
 Gutiérrez [jorgesgutier...@unex.es] Sent: Tuesday, December 10,
 2013 9:48 AM To: r-sig-phylo@r-project.org Subject: [R-sig-phylo] 
 fitting cuadratic model using gnls
 
 Hello everyone,
 
 
 
 I’m trying to fit nonlinear models using the function “gnls” in ape. 
 For example, I run a simple quadratic regression of moult duration 
 (duration) and latitude (latitude_used) as follows:
 
 
 
 mod - duration ~ latitude_used +I(latitude_used^2)
 
 init-list(alpha=1,beta=1)
 
 m3 - gnls(mod, data=data, start=init, correlation = dataCor)
 
 
 
 However, I get the following error message:
 
 
 
 “Error in gnls(mod, data=data, start = init, correlation =
 dataCor) :
 
 step halving factor reduced below minimum in NLS step”
 
 
 
 I’ve also tried with several correlation structures and different 
 starting values for estimation, but without success. I would be very 
 grateful if you can provide any suggestions.
 
 
 
 Thanks in advance,
 
 Jorge S. Gutiérrez
 
 
 [[alternative HTML version deleted]]
 
 
 ___ R-sig-phylo mailing 
 list - R-sig-phylo@r-project.org 
 https://stat.ethz.ch/mailman/listinfo/r-sig-phylo Searchable archive 
 at http://www.mail-archive.com/r-sig-phylo@r-project.org/
 

-BEGIN PGP SIGNATURE-
Version: GnuPG v1.4.11 (GNU/Linux)
Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/

iQEcBAEBAgAGBQJSp2DoAAoJEOCV5YRblxUHn4gH/jYtxFo6aj2DlHAPw2FXufA4
CwtJmIHpWlnD5mhg5YqBgJNldn2XPxqx/qEf0mMzRgd2gAyIQiXd91WI6sqkaVhF
FmYSPzClqx14jZSlF64KE54l/8W4mp81rDyYnh+/SSebltCrX8Cv7NqYqvNDhYQC
M7i/avh0NNv/PZuDuUwbj95doQmdeUyeSlig0sxHaP14C7FxzvxODUBhcTB6yQzN
0wuK+iIaUKqSGm3jmcpY08oy0RJ34QfyH175UezRmSaE/RI+XCohhBbdtU5+WgRX
lBRYmwtCIGAJ+DNmhlnhE3u6ZFdnzhw95gYd1LGL3Sr+JZl1U+VIKfIjhlrTgLQ=
=lLum
-END PGP SIGNATURE-

___
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/