8
The first place to start is the GFF3 specification. This is the official word on what is and is not allowed in a GFF3 file. For example, users can define arbitrary attribute keys, so long as they do not begin with an uppercase letter (these are reserved for "official" use).
But your question doesn't seem to be about what is allowed, but what is commonly ...
5
If I understand the question correctly, you'd like to plot the positions of the matches to you motif along with a gene model that shows the positions of introns and exons for the different transcripts. This can be accomplished fairly easily with ggbio:
library(EnsDb.Mmusculus.v79)
library(ggbio)
library(biomaRt)
library(stringr)
library(dplyr)
# plot gene ...
4
I quite like BED and GFF3 (I don't like GTF/GFF2, though). As text-based formats, I don't think they leave us much room for improvement. Anyway, if you want a new format, here is one. The following is a hybrid between GFF3 and BED. It is a TAB-delimited text-based format with the following fields:
chr, required
start (0-based), required
end, required
strand
...
4
The problem is that GFF, fundamentally, is a relational format: it provides tags that relate individual rows via one-to-many relationships (e.g. gene–exon). This indirectly highlights the second complication: individual rows have different types, and therefore store different attributes in the 9th column.
Over the last few decades (!), we have accumulated a ...
4
Presuming we consider "human readable", "easily parsable", and "quickly queryable" to be objectively good qualities (and if not, I worry for the future):
Text-based: It's absurdly common to want to use grep or awk on annotations. Sure, one could make variants of these that are binary-format aware, but why reinvent the wheel. Of course text files don't ...
4
The three # are used for splitting group of features that belong together, e.g. a transcript and it's exons. Sometimes you see with a blank line instead of the #.
To grep just the header, which has 1 or 2 # from the line start on, you can use extended regular expression:
$ grep -E "^#{1,2}[^#]" Caenorhabditis_elegans.WBcel235.95.gff3
Which means: "grep ...
4
There are some good answers so far, but I don't think any of them fully communicate the significance of the ### directive. The GFF3 specification states:
This directive (three # signs in a row) indicates that all forward references to feature IDs that have been seen to this point have been resolved. After seeing this directive, a program that is ...
3
I stumbled across this issue after getting the same warning message from liftOver. By using the -out option of lgHgGene, the database parameter is ignored but must still be provided. The table parameter should be the second column from your GTF file. Note that you must also specify the -gtf option (otherwise it expects GFF format). In my case, the GTF file ...
3
Usually the home of a program is where the executable is. In order to make the system aware you need to include it as a variable when the terminal is run.
So you need to add it to the file at ~/.bashrc of your computer (if you are using bash as shell [you can check it with echo $SHELL]). You need to add :
export RATT_HOME=path/where/RATT/is/installed/or/...
3
As others have stated, those are just there to separate the entries for easier parsing. They enable you to do nifty tricks like:
$ awk -v RF='###' '/Y74C9A.3/' Caenorhabditis_elegans.WBcel235.95.gff3
I WormBase mRNA 4116 10230 . - . ID=transcript:Y74C9A.3;Parent=gene:WBGene00022277;Name=Y74C9A.3;biotype=protein_coding;transcript_id=...
3
You add 1 to the count of each COG category. You are looking for over represented COG categories so you must count them all. Many genes will be assigned to multiple categories. In fact, a majority of human genes are given multiple GO cellular function annotations. I have worked considerably in this field and one of the surprises I found is how common protein ...
3
I'm one of the EMBLmyGFF3 developers. If the compressed annotation is 396MB as mentioned by @terdon it is not a small annotation... In the worse case we have tested we never waited more than few minutes. We haven't done any particular speed test. It would be valuable if you could tell us how long it has taken in your case (If you didn't kill the process ...
2
That error is telling you that the script you were running was killed by signal 9, also known as SIGKILL, while executing the command _seqret "$@". This is the way the operating system will kill processes that are asking for more memory than it can provide.
Note that removing _secret "$@" is not the solution since that's what's actually doing the job. $@ ...
2
To generate flux balance, one should now the stoichiometric numbers of all the reactions that are taking place, as well as the conditions (that could affect the reactions). I am not sure we have enough knowledge, time and computationally power to do this.
At the moment I'm only aware of a paper simulating the Mycoplasma genitalium cell. But I'm not sure it ...
2
You need to setup ssh keys. All the commands in the script are going to use ssh to run the commands. Note the tips here of how to setup ssh keys:
http://genomewiki.ucsc.edu/index.php/Parasol_job_control_system#SSH_keysssh key setup
2
In A2M format, upper case letters represent matches, lower case letters represent inserts, dashes represent deletions, and dots (or spaces outside the identifier lines) represent gaps aligned to inserts.
So, both "-" and "." are essentially gaps, but assumed to have different origin. This information is supplementary, and most MSA ...
1
that was solved by changing seq level:
newnames <- paste0(c("1","2","3","4", "5","6","7", "8","9","10","11","12","13","14","15","16","17","18","19","20",&...
1
That's a warning. Presumably, a match failed ($& is the default Perl variable which holds the result of a match operation). So something in your input is not exactly what the script was expecting. You can probably safely ignore these messages.
Alternatively, update your question and add a few example lines of input which we can use to reproduce the ...
1
It is possible that the shell being used in the process is different from your standard command line shell. I am not familiar with the tool, and it's hard to say more without knowing your OSX version, but it looks like it is using zsh. I think that zsh only fairly recently became the standard OSX shell. I don't know the tool or the shell very well, but it ...
1
The info is actually in the gene_biotype column of my wormGenes data frame. I could select the rRNA genes using dplyr::filter(wormGenes, gene_biotype == 'rRNA'). Unfortunately, no description is available for the genes, so I will have to get it using the gene IDs.
1
In case an answer would still be relevant after a year and a half, below is a solution with the Biostrings package.
library(Biostrings)
# loading data.table for the fread() function.
# fread() is not only fast but is quite good at guessing how
# the file is formatted
library(data.table)
# readDNAStringSet() function creates a "DNAStringSet" object
...
1
RATT seems to be a very simple wrapper script that unfortunately does not clean up its temporary files properly. Short of modifying the code yourself, your only option is to perform the cleanup manually after running RATT.
The temporary copies are made because (a) the query is not an exact copy: RATT performs some cleaning. And (b) it needs to extract a ...
1
I assume you have the sequence in seq1:
# Set the input
seq <- getSequences("GRIN1")
motif <- "ggcc"
# Convert to lower case
seq2 <- tolower(seqinr::s2c(as.character(seq)))
istarts <- seq(from = 1 , to = length(seq2), by = 1)
# Create words of the same size as the motif
oligos <- seq2[istarts]
for (i in 2:nchar(motif)) {
oligos <- ...
Only top voted, non community-wiki answers of a minimum length are eligible