gmane.comp.lang.r.phylo

length of attribute (names) when calculatingindependent contrasts

Agus Camacho — 2012-05-24T18:58:05

Hi all,

When calculating independent contrasts of a variable, I get the following
message:

  el atributo 'names' [10] debe tener la misma longitud que el vector [9]
The attribute "names" [10] must have the same length than the vector[9]
(this was translated by me)

After checking the help on the command, i made the following

observations,
 > length(Weight)[1] 10> class(Weight)[1] "numeric"> tree$tip.label
[1] "C.leiolepis"      "C.nicterus"       "C.sinebrachiatus"
 [4] "S.catimbau"       "N.ablephara"      "P.erythrocercus"
 [7] "P.tetradactylus"  "V.rubricauda"     "M.maximiliani"
[10] "P.paeminosus"    > class(tree)[1] "phylo"

       0.9020000        0.7155556        0.5370526        0.5887500
     N.ablephara  P.erythrocercus  P.tetradactylus     V.rubricauda
       0.5460000        0.4127273        0.3970000        0.4725294
   M.maximiliani     P.paeminosus
       0.4833333        0.4677692



However, to me, everything seems to be right about the names and the
vector's length. Any hint?

Thanks to all. I am new to these analyses and appreciate much the help
from the mail list.


Agus

problem with fitcontinuos function

Agus Camacho — 2012-05-23T20:16:54

Dear all,

fitcontinuous is giving me a warning generally related to non
correspondence between tree tips and data. However name.check and the tree
tips seems ok.

Anybody has a hint on this?

Thanks in advance.
Agus

Fitting  BM model:> tree
Phylogenetic tree with 11 tips and 10 internal nodes.

Tip labels:
C.leiolepis, C.nicterus, C.sinebrachiatus, S.catimbau, N.ablephara,
P.erythrocercus, ...
Node labels:
1, 2, 3, 4, 5, 6, ...

Rooted; includes branch lengths.

General question about the impact of sampling bias incomparative and diversification analyses

Julien Lorion — 2012-05-23T17:14:21

Dear all,

I have a 5-gene phylogenetic tree, dated, quite well resolved, including 50 species and covering 80~90% of the known diversity... I have also LASER diversification analyses and bayesian reconstruction of ancestral states for 5 characters , both taking into account the impact of sampling biases...

I evaluate now the adaptive significance of 5 characters with ape/caper (using Pagel's test, GEE and pGLS) AND their impact of diversification (using the Yule model with covariates). My results makes perfectly sense in terms of biology but I am worrying about their their relevance. Indeed, for 2 of my characters, I have only 50% of data... In another hand, the distribution of data is well balanced among lineages... and I would say that corresponding patterns of evolution inferred from Bayesian analyses are pretty clear...

I now the Beast team is working on an implementation of a Bayesian model for correlation analysis, that takes into account uncertainties (and thus missing data). But it is not ready yet.

So I have a couple of -likely naive- questions:
1. Is there a way to evaluate the impact of this kind of sampling bias on comparative analyses and diversification analyses ?

2. My guess is that the biases will affect much more directly the Yule model with covariates than GEE, pGLS and Pagel's models... as removing taxa for which I don't have data from the tree will induce an additional bias on diversification rates. Am I right ?

So far, I have removed all missing data in analyses including those 2 characters... My next question may sound crazy and/or desperate and I have a guess on the answer already :-), but...
3. How about assuming that the patterns inferred in my reconstruction of ancestral states are correct (free from unsampled reversal or acquisition that may drastically challenge conclusions) and "completing" real data with most probable states sampled during Bayesian approach along tip branches ? and then compare with results obtained without missing data ?

4 Should I definitely forget about the two characters for which I have only 50% data?

Thanks for any answer/comment/suggestion !

Regards

Julien

Julien Lorion
PhD, Post-doctoral fellow of the Japan Society for the Promotion of Science
Japan Agency for Marine-Earth Science and Technology (JAMSTEC)
Marine Ecosystems Research Department
2-15 Natsushima, Yokosuka 237-0061 Japan
Phone: +81-46-867-9570, Fax: +81-46-867-9525

PIC or PGLS for genome-wide SNP screening

Mattia Prosperi — 2012-05-23T15:05:31

Dear all,

I am working on a data set composed of bacterial genomic sequences (a
few genes) associated to phenotypic values (in-vitro resistance to
antibiotics, a numerical value discretised into a binary class). Of
note, the bacterial isolates were sampled non-uniformly at different
times and locations, thus with a possible sampling bias. The data set
is ~1,000 variables and ~1,000 observations.

I have been applying several methods for developing a model to predict
antibiotic resistance from the single nucleotide polymorphisms (SNP)
extracted from a multiple alignment, applying classical statistical
learning and feature selection methods.
Eventually, I found that a logistic regression with main effects,
where the variables were selected first by a univariable chi-square
screening and then by AIC stepwise, was as good as other more complex
and non-linear methods (such as random forests) by comparing different
loss function (AUROC, specificity, sensitivity) distributions  upon
multiple cross-validation runs. Also, the SNP sets selected by
different approaches were highly similar and consistent across several
bootstrap evaluations.

I found that a few relevant (even after Bonferroni correction) SNP
were located in gene regions that are not supposed to be related with
antibiotic resistance. I thought that this might be a consequence of
neutral mutations that became fixed in the population by chance after
a genetic bottleneck (e.g. antibiotic pressure).
I'd like to understand if such SNP that is associated to antibiotic
resistance (and actually not expected to be) is indeed just a random
mutation of an early isolate that was carrying the true resistance SNP
(in another gene region) and that was selected by the antibiotic
pressure, thus transfering both the true resistance SNP and the
"hitchhicking" ones to the offspring. Unfortunately it is not easy to
cross-tabulate SNP in different genes because not all isolates have
been sequenced the same set of genes.

In order to check for fake/true SNP associated to resistance, I
thought I might use a PIC or PGLS approach (after estimating a
phylogenetic tree from the multiple alignment), in the same settings
as the original analysis, i.e. a model selection approach with both
feature and performance evaluation (well, since the coefficients of
PGLS/OLS are the same, it's just a matter of standard errors and
feature set selection), regressing the resistance class as a dependent
variable and using the SNP as covariates.

Is this a reasonable approach? Does it make sense to set up -for
instance- an AIC stepwise selection within a PGLS modeling?
I know that there is a way to check for phylogenetic signal and
therefore decide if the PGLS approach shall be employed. Is this the
way in which one decides for OLS vs. PGLS?

Last but not least, which is the most appropriate covariance matrix
calculation and PGLS implementation for this input-output set (i.e.
categorical variables, binary class)? The "brunch" function within
caper, or compar.gee within ape?

Thanks and apologies if some of the questions are silly.

M. Prosperi.

Zero value of a significant model parameter duringPGLS

Agus Camacho — 2012-05-23T01:54:25

Dear all,

I just had the following result fitting a Pgls with both, brownian and O-U
correlation structures.
I was wandering what could have made the value and Sd go to zero. Altitude
values were considered continuous but were strongly grouped in two levels.
Might that have to do with the result?

Would appreciate any hint.
There goes the model's summary.
Thanks in advance.


Generalized least squares fit by maximum likelihood
  Model: temperature ~ 1 + altitude
  Data: tp
  AIC BIC logLik
   60  62    -26

Correlation Structure: corMartins
 Formula: ~1
 Parameter estimate(s):
alpha
   21

Coefficients:
            Value Std.Error t-value p-value
(Intercept)    35       2.8    12.7   0.000
altitude        0       0.0    -2.2   0.057

 Correlation:
         (Intr)
altitude -0.93

Standardized residuals:
   Min     Q1    Med     Q3    Max
-1.907 -0.519  0.034  0.734  1.859

Residual standard error: 2.7
Degrees of freedom: 11 total; 9 residual

compute HKY85 distances in R

Juan Antonio Balbuena — 2012-05-22T09:49:55

Hello
I would like to know whether it is possible to compute HKY85 genetic
distances using R. The model is not included in the dist.dna function
of ape and I was wondering whether it is available in another package.
Thank you very much for your attention.

Juan Antonio Balbuena

unrooted tree, spread tips

Walter, Mathias — 2012-05-18T09:42:43

Hi,

I use ape to plot my phylogenetic trees. Often I have to plot unrooted
trees with few clades, but a lot of tips per clade. The tip labels are
hard to read, even if I use lab4ut="axial". Is there any way to spread
the tips? I saw this in some other tree drawing programs.

--
Kind regards,
Mathias

Mentorship program the 2012 iEvoBio conference

Hilmar Lapp — 2012-05-16T18:18:05

The iEvoBio Mentorship Program pairs undergraduate and graduate students with faculty who have extensive experience in evolution, systematics, biodiversity, mathematics and/or software development to enhance their experience at iEvoBio 2012. Mentors will guide participants in collaborative and networking opportunities throughout the two-day conference.

If you would like to gain the wisdom of your colleagues' experiences by participating in this program, please email the iEvoBio 2012 organizing committee at committee-Phy/wAXTOfxAfugRpC6u6w< at >public.gmane.org by June 1, indicating if you have any activities you would like to be included in the mentorship program, or a particular field of interest. You can find more information about the meeting and registration at http://ievobio.org/.

iEvoBio 2012 is sponsored by the US National Evolutionary Synthesis Center (NESCent) and by Biomatters Ltd., in partnership with the Society for the Study of Evolution (SSE) and the Systematic Biologists (SSB).

The iEvoBio 2012 Organizing Committee:
Hilmar Lapp, US National Evolutionary Synthesis Center (chair)
Robert Beiko, Dalhousie University
Nico Cellinese, University of Florida
Robert Guralnick, University of Colorado at Boulder
Rebecca Kao, Denver Botanic Gardens
Ellinor Michel, Natural History Museum, London
Nadia Talent, Royal Ontario Museum
Andrea Thomer, University of Illinois at Urbana-Champaign

Best way to test correlation between discrete andcontinuous variables ?

Julien Lorion — 2012-05-16T07:53:38

Dear all, 

I am working on the evolution of deep-sea symbiotic mussels... I have got a tree and 5 characters: habitat (hydrothermal vents, cold seep and organic substrate), presence/absence of methanotrophic symbionts, presence/absence of sulfoxydizing symbionts, symbiont location (extra VS intracellular) and body length... 

So that's 1 continuous and 4 discrete binary variables (actually, I assumed vent and seeps are very similar... so no need to take into account the 3 states) 

At first I tested various hypotheses about correlation between my discrete characters... I chose the easy way: I remembered my master lectures and used basic Pagel's correlations. If you think that any new tool performs better, I'd be happy to hear it. 

For now, my main concern is that I wanna test the impact of two binary variables (habitat and symbiont location) on the body length... In fact, that looks likes ANOVAs from which I wanna remove the phylogenetic bias... 

The point is that I don't know how to do that in practice. May someone have some advices ?

Thanks by advance

Best regards

Julien

exporting nexus tree in APE

Lucinda Kirkpatrick — 2012-05-15T00:38:51

Dear all,

I am attempting to export a tree I have pruned from 916 taxa to the 53 taxa
I am interested. I used this script:

overlap <- name.check(tree, d)
tree <- drop.tip(tree, overlap$Tree.not.data)
d2 <- d[!(rownames(d) %in% overlap$Data.not.tree),]


plot(tree)
name.check(tree,d2)

That was successful but now I want to export the reduced tree as a nexus
file I get the following error:

 write.nexus(tree, file="./Results/Tree.nex")
Read 921 items
Error in 1:(start - 1) : argument of length 0

I need the reduced tree in a nexus format for another analysis - how is the
best way to do that?  I am very new to R so please make any explanation
very simple!

Thanks, Luci

[[alternative HTML version deleted]]

NaN returned with dist.dna {ape}

Augusto Ribas — 2012-05-14T21:25:19

Hello.
I'm reading" Analysis of Phylogenetics and Evolution with R" and tryed to
make a tree with some frogs species repeating the Sylvia Case study example.
I had no problem in search, download and align many 16s sequences for
several frog species.
I used a web service to align the sequences.
The results are here:
http://www.ebi.ac.uk/Tools/services/web_clustalw2/toolresult.ebi?tool=clustalw2&jobId=clustalw2-I20120514-220954-0303-43791553-oy

Then i read this align sequences and try to calc the distance, but some NaN
are returned, not allowing me to continue.

#reading data
sapinhos.seq.ali<-read.dna("
http://www.ebi.ac.uk/Tools/services/rest/clustalw2/result/clustalw2-I20120514-220954-0303-43791553-oy/aln-phylip
")
sapinhos.seq.ali

#caling distances, some NaN return for every distance
sapinho.K80<- dist.dna(sapinhos.seq.ali, pairwise.deletion = T)
sapinho.F84<- dist.dna(sapinhos.seq.ali, model = "F84",pairwise.deletion =
TRUE)
sapinho.TN93<- dist.dna(sapinhos.seq.ali, model = "TN93",pairwise.deletion
= TRUE)
sapinho.GG95<- dist.dna(sapinhos.seq.ali, model = "GG95",pairwise.deletion
= TRUE)

round(cor(cbind(sapinho.K80,sapinho.F84,sapinho.TN93,sapinho.GG95)),3)


Nothing worked from here. I don't understand what is wrong.
Can someone lend me a hand.

PS. According to the site, the data is saved for 7 days, so if this e-mail
past 7 day we probably wont be able to read the data direct from the site,
but i have every step saved here if anything else can help
in the solution, just tryed to be short.

Distance Matrix for morphological characters, including "?", can it still be euclidean?

Robert Schenck — 2012-05-14T03:11:26

I've been trying to create a distance matrix from a dataset of
morphological characters, where the matrix contain some characters as
"?". Is it possible to get a specifically euclidean matrix from that?
When I try to make the distance matrix, I get a what looks to me like
a distance matrix,, but it warns that NAs were created by coercion.
So I am doing: "dist(x, method="euclidean")" but "is.euclid(x)"
responds that it's not.

pGLS with populations from diferent localities

Agus Camacho — 2012-05-13T20:30:33

Dear all,

I am trying to compare environmental vs phylogenetic effects over one trait
that has evolved in a single branch. For that, I have 10 closely species
that have been collected over 7 localities. The number of localities in
which each species has been collected varies from 1 to 3. So some species
are endemic and others are more widespread among localities.

I am wondering how should I do to relate the data set and the tree: Is it a
nonsense to repeat the name of the widespread species in the tree?
Would appreciate any hint on this.
Thanks in advance.
Agus

Calculation of NRI/NTI for large species matrix

micro buggy — 2012-05-13T16:57:21

Dear folks,

I tried to calculate NRI/NTI for large community matrix (species
larger than 20,000, for instance) using Picante. My current PC memory
size (12 Gb) can not hold this kind of analyses. Does someone know how
to do with this issue? Thanks.

Si

nonparametric PGLS

Pascal Title — 2012-05-11T03:46:32

Hi everyone

I have character data for species on a phylogeny, but simple
transformations like log-transformations or square-root transformations are
not proving to be sufficient to get the data to be normal. If this were a
non-phylogenetic test, I would resort to something like a Spearman
correlation, which is non-parametric.
But with phylogenetic generalized least squares, is there any method that
is nonparametric?

Thanks!


-Pascal Title

[[alternative HTML version deleted]]

OUwie v.1.21

Jeremy Beaulieu — 2012-05-09T15:41:25

Hi Everyone~

We have just released a new version of OUwie (v.1.21), which is now available on CRAN (and R-Forge). OUwie provides an implementation of a new set of Ornstein-Uhlenbeck-based Hansen models that allow the strength of selection, optimal trait value, and/or stochastic motion parameters to vary across a tree. Despite the small change in version number since the initial release, this has a few major new features: 1) Trees with mapped character reconstructions saved in SIMMAP format can now be loaded into OUwie and the different character states painted with different models; 2) We now allow users to define the bounds on the parameter search; 3) Users can easily specify a new selective regime for a particular clade of interest without having to paint a tree manually; and 4) Overall speed has bee
 n increased nearly three-fold. 

There will be continued improvements to OUwie, and over the coming months we plan to release several new functions and capabilities. Rather than continue to use R-sig-phylo for major announcements, we have set up two OUwie mailing lists: http://groups.google.com/group/ouwie-announce for announcements of new features, and http://groups.google.com/group/ouwie-discuss as a discussion forum for OUwie problems, suggestions, and other issues. Thank you for the feedback we have received so far -- it has been quite helpful.

Sincerely,
Jeremy Beaulieu
Brian O'Meara

_________________________________
Jeremy Beaulieu
Graduate Student
Dept. Ecology & Evolutionary Biology
Yale University
www.jeremybeaulieu.org

Problems useing query() function of package "seqinr"

Augusto Ribas — 2012-05-09T14:00:02

Hello, i'm new in this list, and new in genetics field too.
I'm reading a book "Analysis of Phylogenetics and Evolution with R" and in
the book i saw about the function query() of the package seqinr to search
databases about specific sequences of  our interest.

#A got a list of some frog species:

lista.spp<-c("Rhinella rubescens", "Rhinella fernandezae", "Elosia
rustica",
"Crossodactylus gaudichaudii", "Leptodactylus pustulatus", "Leptodactylus
syphax",
"Dendropsophus cachimbo", "Leptodactylus marmoratus", "Physalaemus
soaresi",
"Hylodes nasus", "Scinax fuscomarginatus", "Leptodactylus petersii",
"Chiasmocleis capixaba", "Ceratophrys aurita", "Pleurodema diplolister",
"Cystignatus gigas", "Scinax trilineatus", "Elosis nasus", "Dryadophis
bifossatus"
)

#And stared looking for 16s RIBOSOMAL RNA.
#My problem is that for some species it worked great.
#For example:

[1] 1
[1] "GU907196.RR1"

#1 sequence and the name to get the sequence later

[1] 2
[1] "JF789923" "JF789924"

#2 sequences and names to get them later.

#The problems is that for some species like

R stops, i thought it would return too many information, but even using the
argument virtual=T it stops the same way, and i have to turn R off.
I looked in the genbank site and this species should return no sequence but
i don't know what is wrong.
Also, in my understanding, "Leptodactylus syphax" should return 2
sequences, looking in genbank site, but it also stops.
But from the 19 species, only some stops, the others work, i cant figure
out what is happening.

I use windows 7 OS, a residential internet, R 2.15.0 and Tinn-R as a R-gui

Also, in the university the port query() function use is blocked, only
ports like 80 are open in the firewall, are there other functions like this
in R?
In the task view, biocondutor project was cited but i still searching there.

Could someone enlighten me?

Best Wishes
Augusto Ribas

estimate ancestral states of specific internal nodes

Annemarie Verkerk — 2012-05-08T23:44:25

Dear r-sig-phylo users,

I have a question regarding the retrieval of ancestral states of 
internal nodes in a sample of trees. I am doing ancestral state 
estimations on a large sample of trees (using ace() in ape), but a few 
groups of languages are present in all my trees. I would like to know 
how to easily retrieve the ancestral state of the internal node leading 
to such a group of languages, given that I know the numbers of the edges 
leading to those language tips, but the number of the internal edge 
leading to those languages (and thus defining the grouping / clade) will 
not be the same in all of the trees.

This would be very easy to do in BayesTraits using 'addnode', but I 
really want to use the 'ML' version of ace(), which is not included in 
BayesTraits.

Does anyone have an idea? (Sorry if this has been asked before, I could 
not find anything in the archives.)

Many thanks,
Annemarie

pic.ortho data entry

Benstead, Jonathan — 2012-05-08T16:55:05

Dear r-sig-phylo members - Some colleagues and I are trying to apply the pic.ortho function in ape, but we're struggling with how to enter the data for the analysis. We have 72 data points from 37 species. We know that the input data have to be in the form of a list of numeric vectors, one vector for each species, but we're unsure about how these vectors get matched up with the information from the tree. There doesn't seem to be much on the web, and the example in the manual uses randomly generated data. Any help would be greatly appreciated. Jon.

compar.gee: how to extract SE, t, and P from output?

Nina Hobbhahn — 2012-05-08T11:55:24

Dear fellow list users,

I would like to extract the values for SE, t, and Pr from the compar.gee output, but seem to unable to do so. The command lines

GEE<-Ââcompar.gee(trait~reward-Ââ1, phy=phy2, data=DF.Disa, family=gaussian) 

GEE

return the following output:

Beginning Cgee S-Ââfunction, < at >(#) geeformula.q 4.13 98/01/27 running glm to get initial regression estimate

rewardN rewardY
2.689142 3.421509
Beginning Cgee S-Ââfunction, < at >(#) geeformula.q 4.13 98/01/27 running glm to get initial regression estimate

rewardN rewardY
2.689142 3.421509
Call: compar.gee(formula = trait ~ reward -Ââ 1, data = DF.Disa,

family = gaussian, phy = phy2) Number of observations: 31

Model:
Link: identity

Variance to Mean Relation: gaussian


QIC: 165.2011

Summary of Residuals:
Min 1Q Median 3Q Max

3.1462585 0.2560335 1.5256031 2.2012787 4.6689764

Coefficients:
Estimate S.E. t Pr(T > |t|)

rewardN 1.959375 1.091801 1.794626 0.1016430 rewardY 1.908384 1.117150 1.708261 0.1170666

Estimated Scale Parameter: 4.862036 "Phylogenetic" df (dfP): 12.45332

When typing 

GEE$coefficients

I get


rewardN rewardY

1.959375 1.908384


but I'm unable to access SE, t, and Pr(T > |t|). Please can any of you help? Thank you very much,


Nina




Dr. Nina Hobbhahn
Post-doctoral fellow
Lab of Prof. S. D. Johnson
School of Life Sciences
University of KwaZulu-Natal
Private Bag X01
Scottsville, Pietermaritzburg, 3201
South Africa
[[alternative HTML version deleted]]