Data Investigation - Mysterious Case of Excessive Exons in CEABIGR
CEABIGR
Crassostrea viriginica
Eastern oyster
2023
Author
Sam White
Published
December 6, 2023
After Steven ran through some prelimiary exon expression analyses (GitHub Isse), he noticed that there were some genes with a seemingly excessive number of exons (like over 800!). Seeing that, I’ve decided to do some investigation to see if these gene feature annotations stem from the source (i.e. from NCBI GFF) or were introduced accidentally through some subsequent data analysis/manipulations.
First, let’s look at a file which appears to provide mappings of GeneIDs and exons:
Well, this shows us that this file indeed has many, many genes with thousands of exons! This is highly dubious, so let’s glance at the NCBI GFF to see if this is “real” or not.
So, the exon counts in genome-features/C_virginica-3.0_Gnomon_exon-geneID.bed are indeed artificial!
Now that I’ve investigated that, I probably should’ve just started with the input files that Steven was using in the first place. However, this certainly identifies a problematic file that should likely be removed from the repo…
Let’s look at the data Steven was working with and see what we can figure out.
Count colums in output/19-exon-expression/S12M-exon_expression.tab
Due to the layout of this file (a column for each possible number of exons across all genes), the head command will produce an unwieldy output. Instead, we’ll just count the number of columns (fields) to get an idea of what we’re dealing with.
'data.frame': 38838 obs. of 463 variables:
$ gene_name: chr "gene-LOC111116054" "gene-LOC111126949" "gene-LOC111110729" "gene-LOC111112434" ...
$ exon_1 : int 130 22 282 24 69 416 193 23 23 3 ...
$ exon_2 : int 103 180 63 40 39 NA 106 173 173 1 ...
$ exon_3 : int 9 35 55 123 86 NA NA 3 NA 2 ...
$ exon_4 : int NA 110 11 NA 21 NA NA 1 NA 2 ...
$ exon_5 : int NA 47 NA NA 85 NA NA 2 NA NA ...
$ exon_6 : int NA NA NA NA NA NA NA 2 NA NA ...
$ exon_7 : int NA NA NA NA NA NA NA 75 NA NA ...
$ exon_8 : int NA NA NA NA NA NA NA 214 NA NA ...
$ exon_9 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_10 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_11 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_12 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_13 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_14 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_15 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_16 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_17 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_18 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_19 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_20 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_21 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_22 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_23 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_24 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_25 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_26 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_27 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_28 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_29 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_30 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_31 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_32 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_33 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_34 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_35 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_36 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_37 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_38 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_39 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_40 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_41 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_42 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_43 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_44 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_45 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_46 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_47 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_48 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_49 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_50 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_51 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_52 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_53 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_54 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_55 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_56 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_57 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_58 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_59 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_60 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_61 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_62 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_63 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_64 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_65 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_66 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_67 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_68 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_69 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_70 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_71 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_72 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_73 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_74 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_75 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_76 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_77 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_78 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_79 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_80 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_81 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_82 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_83 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_84 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_85 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_86 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_87 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_88 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_89 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_90 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_91 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_92 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_93 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_94 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_95 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_96 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_97 : int NA NA NA NA NA NA NA NA NA NA ...
$ exon_98 : int NA NA NA NA NA NA NA NA NA NA ...
[list output truncated]
result <- S12M.df %>%filter(!is.na(exon_462) &grepl("[0-9]", exon_462))print(result$gene_name)
[1] "gene-LOC111119012"
selected_gene_names <- result$gene_name
Okay, we’ve identified gene-LOC111119012 as the gene which contains 462. Let’s see if this matches the original NCBI GFF.
Compare number of exons in gene-LOC111119012 between S12M-exon_expression.tab and GCF_002022765.2_C_virginica-3.0_genomic.gff
Okay, they do not match, which points to an issue with S12M-exon_expression.tab! I generated that file in 19-exon-expression.Rmd using the Ballgown exon file (e_data.ctab) file as a source. Let’s check that first before I re-examine 19-exon-expression.Rmd.
Compare total exons in e_data.ctab and GCF_002022765.2_C_virginica-3.0_genomic.gff
So, we’re seeing that S12M e_data.ctab has ~50% of the total number of exons in the Crassostrea virginica (Eastern oyster) genome (as annotated by NCBI). This makes sense to me, as the data we’re using is RNAseq. So, not all genes are going to be expressed and, thus, we won’t have data for each gene/exon in the genome.
Regardless, based off of the difference in exon counts of gene-LOC111119012 in GCF_002022765.2_C_virginica-3.0_genomic.gff and S12M-exon_expression.tab, I will revisit 19-exon-expression.Rmd to see where things might have gone wrong.