cgeb2001's DokuWiki!

User Tools

Site Tools

Trace:
gene_prediction_curation_with_igv

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision		Previous revision
gene_prediction_curation_with_igv [2023/10/18 09:59] – [Validating curated gene models] 134.190.232.191gene_prediction_curation_with_igv [2025/08/07 16:01] (current) – 134.190.145.228
Line 26: Line 26:
 Note that here two distinct RNAseq tracks are loaded. One for all the reads stemming from transcripts transcribed from the positive strand of the genome, and one for all those transcribed from the negative strand of the genome. This kind of strand-specificity information is only available when a RNAseq library protocol was used to retains this information, for example the dUTP method. I believe nowadays this is the default so it should be available to you. Note that here two distinct RNAseq tracks are loaded. One for all the reads stemming from transcripts transcribed from the positive strand of the genome, and one for all those transcribed from the negative strand of the genome. This kind of strand-specificity information is only available when a RNAseq library protocol was used to retains this information, for example the dUTP method. I believe nowadays this is the default so it should be available to you.
  
-If you used HISAT2 to map your RNAseq reads to your genome, all "positive" reads will have the tag ''XS:A:+'', while all "negative" reads will have the ''XS:A:-'' tag. One can separate the two types of reads from the HISAT2-generated BAM file using the following command:+If you used HISAT2 with ''--rna-strandness RF'' to map your RNAseq reads to your genome, all "positive" reads will have the tag ''XS:A:+'', while all "negative" reads will have the ''XS:A:-'' tag. One can separate the two types of reads from the HISAT2-generated BAM file using the following command:
  
 <code> <code>
Line 35: Line 35:
 samtools index <rnaseq.neg.bam> <rnaseq.neg.bam.bai> samtools index <rnaseq.neg.bam> <rnaseq.neg.bam.bai>
 </code> </code>
 +
 +==== GC content ====
 +
 +It may be nice to have a track in IGV that shows you per-window-size %GC. Here I generate a so-called ''bedGraph'' formatted file that can be loaded into IGV
 +
 +<code>
 +cat Ergobibamus_cyprinoides_CL.scaffolds.fa \
 +    | seqkit sliding -W 50 -s 25 \
 +    | seqkit fx2tab -n -g \
 +    | sed -e "s/_sliding\:/\t/" -e "s/\-/\t/" -e "s/\s+/\t/g" \
 +    > ergo_gc_content.bedgraph
 +</code>
 +
 +You can play with the window size ''-W'' and overlap ''-s'' parameters to fine tune it to your own liking. You can either load in this track via ''File -> Load From File ..'' or through the JSON file (see below)
  
 ==== Loading into IGV ==== ==== Loading into IGV ====
Line 156: Line 170:
             "url": "rnaseq_vs_masked_ergo_cyp_genome.sort.neg_transcripts.bam",             "url": "rnaseq_vs_masked_ergo_cyp_genome.sort.neg_transcripts.bam",
             "indexURL": "rnaseq_vs_masked_ergo_cyp_genome.sort.neg_transcripts.bam.bai"             "indexURL": "rnaseq_vs_masked_ergo_cyp_genome.sort.neg_transcripts.bam.bai"
 +        },
 +        {
 +            "name": "GC content",
 +            "type": "wig",
 +            "format": "bedGraph",
 +            "url": "ergo_gc_content.bedgraph",
 +            "graphType": "bar",
 +            "min": 20,
 +            "max": 80,
 +            "color": "#00dad7"
         }         }
 +
     ]     ]
  
Line 177: Line 202:
 "name": Desired display name for your track "name": Desired display name for your track
 "type": Type of track. Could be "alignment", "annotation" and perhaps other types "type": Type of track. Could be "alignment", "annotation" and perhaps other types
-"format" : "bam" or "gff3"+"format" : "bam" or "gff3" or "wig"
 "url" : filename of the file you wish to load here "url" : filename of the file you wish to load here
 "fastaURL" : filename of the genome FASTA file "fastaURL" : filename of the genome FASTA file
Line 204: Line 229:
  
 For some frustratingly unknow reason, IGV only displays the three forward frames, OR the three reverse frames. It is, as far as I know, not possible to display all six frames simultaneously. To swap between reverse and forward frames, click on small arrow of the //Sequence// track. To see the sequence track, you need to be sufficiently zoomed in. For some frustratingly unknow reason, IGV only displays the three forward frames, OR the three reverse frames. It is, as far as I know, not possible to display all six frames simultaneously. To swap between reverse and forward frames, click on small arrow of the //Sequence// track. To see the sequence track, you need to be sufficiently zoomed in.
 +
 ==== Curating gene models ==== ==== Curating gene models ====
  
Line 330: Line 356:
   * Does the translate CDS start with an M and end with a * (STOP)?   * Does the translate CDS start with an M and end with a * (STOP)?
   * Are there any premature * (STOP) in the CDS?   * Are there any premature * (STOP) in the CDS?
 +  * If Blastocystis genome, and there is no STOP, does it end with a TGTTTGTT motif?
  
 <code> <code>
Line 345: Line 372:
 </code> </code>
  
-NOTE: Currently, the checks performed in the script assume that you haven't annotated any UTRs in your GFF3. If you have UTR (''three_prime_UTR'' or ''five_prime_UTR'') features, for example the start coordinate of the first CDS will not match the start coordinate of the gene. It will thus throw an error, even though your GFF3 is perfectly fine. Something I'll have to fix in the future.+<del>NOTE: Currently, the checks performed in the script assume that you haven't annotated any UTRs in your GFF3. If you have UTR (''three_prime_UTR'' or ''five_prime_UTR'') features, for example the start coordinate of the first CDS will not match the start coordinate of the gene. It will thus throw an error, even though your GFF3 is perfectly fine. Something I'll have to fix in the future.</del> 
 + 
 +UPDATE (Jan 2025): It should now account for UTRs 
 + 
 +==== Re-streamlining gene names ==== 
 + 
 +During your curation, you may be forced to make new gene names that stray from the original format. For example, if you found that gene ''FUN_000012'' was actually two genes that were wrongly fused, you may have split them up into ''FUN_000012-1'' and ''FUN_000012-2''. This is perfectly fine, but may cause problems if you want to feed the curated GFF3 as input to for example ''funannotate annotate'' or ''funannotate update''.  
 + 
 +To re-streamline gene names, you can make use of the ''funannotate util gff-rename'' utility: 
 + 
 +<code> 
 +funannotate util gff-rename \ 
 +    --gff3 Ergobibamus_curated.gff3 \ 
 +    --fasta Ergobibamus_contigs.fasta \ 
 +    --locus_tag PYV62 \ 
 +    --out Ergobibamus_curated_renamed.gff3 
 +</code>
  
-TODO: also make it so that it collects all errors are reports them at the end, rather than quitting at the first error. Make use of assert statements?+Note that I have also added a new locus tag here, ''PYV62''. This was the locus tag assigned to me by the NCBI genome submission portal.
gene_prediction_curation_with_igv.1697633973.txt.gz · Last modified: by 134.190.232.191