[Seqinr-forum] Query regarding running "Online Synonymous Codon Usage Analyses with the ade4 and seqinR packages"
sourav roy choudhury
src.bio at gmail.com
Thu Nov 5 19:43:23 CET 2009
Dear Dr. Lobry and Colleagues,
The code is running fine and the codon count table is generated.
But the problem seems to be that it only produces codon table for the first
virus on the list, NC_000898 and not for all the available viruses. Another
concern is that for the virus NC_000898 have 104 CDS but it seems to be only
6 CDS is selected here, is it so? Additionally I do not want to include "Met",
"Trp" and "Stop" Codons for the correspondence analysis.
A minor addition, I want to select only those CDS with length >150 bp or >50
Codons.
I really appreciate your help. Thank you very much
Sincerely yours
Sourav Roy Choudhury
SOURAV ROY CHOUDHURY.
Dept of Neuroscience
University of Calcutta
West Bengal, INDIA
On Thu, Nov 5, 2009 at 1:48 AM, Jean lobry <lobry at biomserv.univ-lyon1.fr>wrote:
> The code is running fine and I can reproduce the output as
>> same as yours. It seems that the data is accessible.
>>
>
> Perfect.
>
>
> Yes, these are absolutely those genomes that I intend to
>> work with. I want to consider only all the CDS (or Open
>> Reading Frame) for the statistical observation.
>>
>
> OK, let's make one more baby step: can you compute the table
> of codon counts with the following code?
>
> ######
> library(seqinr)
> choosebank("refseqViruses")
> #
> # Now we want all to extract all complete CDS from
> # these genomes:
> #
> query("herpcds", "TID=548681 AND T=CDS AND NO K=PARTIAL")
> n <- herpcds$nelem
> n # there are 4887 CDS
> mnemo.cds <- getName(herpcds)
> head(mnemo.cds) # just to see what they look like
> #
> # We pre-allocate a table for codon counts:
> #
> tabcod <- matrix(0, nrow = n, ncol = 64)
> rownames(tabcod) <- mnemo.cds # cds names
> colnames(tabcod) <- words() # codon names
> head(tabcod) # this empty for now
> #
> # Now we loop over all CDS to compute codon usage (this
> # may take a while)
> #
> pb <- txtProgressBar(min = 1, max = n, initial = 1, style = 3)
> for(i in 1:n){
> setTxtProgressBar(pb, i) # to show progression
> tabcod[i, ] <- uco(getSequence(herpcds$req[[i]]))
> }
> close(pb)
> head(tabcod) # shouldn't be empty now
> #
> # We save this table on disk for further usage
> #
> save(tabcod, file = "tabcod.RData")
> ######
>
>
> Best,
> --
> Jean R. Lobry (lobry at biomserv.univ-lyon1.fr)
> Laboratoire BBE-CNRS-UMR-5558, Univ. C. Bernard - LYON I,
> 43 Bd 11/11/1918, F-69622 VILLEURBANNE CEDEX, FRANCE
> allo : +33 472 43 27 56 fax : +33 472 43 13 88
> http://pbil.univ-lyon1.fr/members/lobry/
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.r-forge.r-project.org/pipermail/seqinr-forum/attachments/20091106/70ab59cb/attachment.htm
More information about the Seqinr-forum
mailing list