Re: [R-sig-phylo] Extract singleton haplotypes with pegas::haplotype()

2021-04-26 Thread Emmanuel Paradis
Hi Jarrett,

- Le 27 Avr 21, à 0:09, Jarrett Phillips phillipsjarre...@gmail.com a écrit 
:
> Hi All,
> 
> I'm looking for a way to extract singleton haplotypes using the haplotype
> function in Pegas.
> 
> I use the below function to accomplish this:
> 
>extract.singleton.haps <- function(freq) {
>  ind <- which(freq == 1) # index of singleton haplotypes
>  single <- seqs[ind, ] # extracted haplotypes (only singletons)
> }

The object 'seqs' is neither created inside the function nor passed as 
argument, so R will look for it using the infamous lexical scoping rule, 
meaning that you're not sure it'll do what you want it to do.

> I have tried the following:
> 
>library(pegas)
>seqs <- read.dna(file = file.choose(), format = "fasta") # select FASTA
> file
>freq <- haplotype(seqs) # haplotype distribution
>freq <- lengths(attr(freq, "index") # extract haplotypes

In the line above, you overwrite 'freq' created in the previous line, so maybe 
you can change these two for:

haps <- haplotype(seqs) # haplotype sequences and distribution
freq <- lengths(attr(haps, "index") # extract frequencies

>singleton.seqs <- extract.singleton.haps(freq)

See my comment above: 'seqs' exists in your global environment so the function 
singleton.seqs will use the original sequences (because of lexical scoping), 
not the haplotype ones. You can do:

singleton.seqs <- haps[freq == 1, ]

or:

singleton.seqs <- seqs[attr(haps, "index")[freq == 1, ], ]

I'm not sure that writing a function is really needed here: probably you want 
to examine all values in 'freq'.

HTH

Emmanuel

>write.dna(singleton.seqs, file = paste0(taxon, "_singletons.fas"),
> format = "fasta")
> 
> 
> However, upon viewing the sequences in MEGA, they appear to all be
> identical, instead of unique.
> 
> Any ideas on what is happening here?
> 
> Thanks.
> 
> Cheers,
> 
> Jarrett
> 
>   [[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/


[R-sig-phylo] Extract singleton haplotypes with pegas::haplotype()

2021-04-26 Thread Jarrett Phillips
Hi All,

I'm looking for a way to extract singleton haplotypes using the haplotype
function in Pegas.

I use the below function to accomplish this:

extract.singleton.haps <- function(freq) {
  ind <- which(freq == 1) # index of singleton haplotypes
  single <- seqs[ind, ] # extracted haplotypes (only singletons)
}

I have tried the following:

library(pegas)
seqs <- read.dna(file = file.choose(), format = "fasta") # select FASTA
file
freq <- haplotype(seqs) # haplotype distribution
freq <- lengths(attr(freq, "index") # extract haplotypes
singleton.seqs <- extract.singleton.haps(freq)
write.dna(singleton.seqs, file = paste0(taxon, "_singletons.fas"),
format = "fasta")


However, upon viewing the sequences in MEGA, they appear to all be
identical, instead of unique.

Any ideas on what is happening here?

Thanks.

Cheers,

Jarrett

[[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] Are geiger:dtt p-values one-tailed?

2021-04-26 Thread Roberto Marquez
Hello everyone. I have a quick question about the p-values calculated by the 
dtt function in geiger ($MDIpVal in the function output): Does the p-value 
correspond to p[MDI>=0] (ie. one-tailed) or p[abs(MDI)>=0] (two-tailed)? 
Looking at the code for getMDIp() it seems like the p-value is one-tailed, but 
I wanted to confirm.

Thanks!

Roberto

Roberto M�rquez
Postdoctoral Scholar
Department of Ecology and Evolution
University of Chicago

[[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] ape: error in pcoa

2021-04-26 Thread Juan Antonio Balbuena
Many thanks Emmanuel for the detailed instructions and sorry for the 
messed up formatting.

Best

Juan


El 24/04/2021 a las 9:03, Emmanuel Paradis escribió:
> Hi Juan,
>
> There is a bug in pcoa() when only one axis is selected by the function. You 
> can fix it by editing the function with fix(pcoa) and modify line 38 (or 
> close depending on how the editor arranges the lines of code):
>
> vectors <- sweep(D.eig$vectors[, 1:k], 2, sqrt(eig[1:k]), FUN = "*")
>
> insert ', drop = FALSE' where D.eig$vectoes is indexed so this line becomes:
>
> vectors <- sweep(D.eig$vectors[, 1:k, drop = FALSE], 2, sqrt(eig[1:k]), FUN = 
> "*")
>
> save and close. There are two other similar lines but they are not used in 
> your example. I've fixed that in ape 5.5 (which should be on CRAN soon).
>
> As a side note, it took me some time to figure out the shape of the data you 
> posted. Mail applications sometimes format lines of text in different ways, 
> so it's better to run a command like this:
>
> cat(deparse(DH), sep = "\n")
>
> and paste the output in your message. This can then be input with: DH <- 
> [[the output from above]] with the attributes correctly set.
>
> Best,
>
> Emmanuel
>
> - Le 23 Avr 21, à 21:24, Juan Antonio Balbuena j.a.balbu...@uv.es a écrit 
> :
>
>> Hellow,
>>
>> I am running a set of simulations with distance matrices used as input
>> to pcoa (package ape).
>>
>> In one instance the matrix (DH) is
>>
>> H18 H20 H21 H18 0.000 0.3127452190 0.3127452190 H20 0.3127452
>> 0.00 0.0001625185 H21 0.3127452 0.0001625185 0.00
>>
>> when I run pcoa I get this
>>
>> H.PCo <- pcoa(DH, correction="none")$vectors Error in array(STATS,
>> dims[perm]) : 'dims' cannot be of length 0
>>
>> However, if I multiply DH*10 pcoa runs just fine:
>>
>> H.PCo <- pcoa(DH*10, correction="none")$vectors
>>
>> A word of explanation will be very much appreciated.
>>
>> Thank you for your attention
>>   
>> Juan A. Balbuena
>>
>> --
>>
>> Dr. Juan A. Balbuena
>> Cavanilles Institute of Biodiversity and Evolutionary Biology
>> Symbiont Ecology and Evolution Lab
>> University of Valencia http://www.uv.es/~balbuena
>> 
>> P.O. Box 22085 http://www.uv.es/cophylpaco 
>> 46071 Valencia, Spain
>> e-mail: j.a.balbu...@uv.es tel. +34 963 543
>> 658    fax +34 963 543 733
>> 
>> *NOTE!*For shipments by EXPRESS COURIER use the following street address:
>> C/ Catedrático José Beltrán 2, 46980 Paterna (Valencia), Spain.
>> 
>>
>>
>>  [[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/
-- 

Dr. Juan A. Balbuena
Cavanilles Institute of Biodiversity and Evolutionary Biology
Symbiont Ecology and Evolution Lab
University of Valencia http://www.uv.es/~balbuena 

P.O. Box 22085 http://www.uv.es/cophylpaco 
46071 Valencia, Spain
e-mail: j.a.balbu...@uv.es tel. +34 963 543 
658    fax +34 963 543 733

*NOTE!*For shipments by EXPRESS COURIER use the following street address:
C/ Catedrático José Beltrán 2, 46980 Paterna (Valencia), Spain.



[[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/