I'll try to tell the histoy short

I have stranded (dUTP) RNA-Seq paired data obtained with Illumina (2x150) from red and green fruits. A total of three biological replicates for each

I ran kallisto quant (b 100, -rf option), using these reads and the latest annotated fasta CDS file which contains some isoforms from many of the genes. The CDS file I indexed can be found here FASTA CDS File. Over 85% of the paired reads mapped to the fasta file

Then, the corresponding h5 abundance files were imported into DESeq2 using tximport, and a DE analysis was performed. That was pretty easy as only two conditions were compared

I was certainly happy to see that in a particular gene of interest, a significant Fold Change value was obtained (green vs red). I am talking of a particular gene that have 5 different isoforms (t1 to t5) described in the fasta file. CDS corresponding to the t1, t2, t4 and t5 isoforms had a fold change ranging from 4 to 6, but I was certainly excited to see that t3 has a Fold Change over 256. This change could be very meaningful for me from a biological point of view.

My first expectation is that the t3 isoform should have some components (i.e. a whole new exon being included in the transcript or the like) that should be missing in the rest of the other isoforms, and that a regulated and specific expression of this t3 isoform took place in the red fruits only. Reads that are mapped to this new component could explain the observed Fold CHange

To my surprise, after running a multiple alignment with Clustal Omega, I saw that t1and t3 are almost identical, being t3 only a little bit shorter in the 3' side (30 to 40 bases shorter). The rest of the sequence were the same

Thus, I have not a reasonable way to explain these results.

Code should be placed in three backticks as shown below

# include your problematic code here with any corresponding output 
# please also include the results of running the following in an R session 

sessionInfo( )

From my lab, we have a method Swish which is a nonparametric method based on SAMseq designed to help with isoform level analysis. When isoforms have similar sequence, it's important to keep track of the quantification uncertainty. We found Swish performed well in these cases compared to parametric models. You just need to have bootstrap or Gibbs posterior samples from the upstream quantification method. Let me know if you need further pointers, Swish is in the fishpond Bioc package:

https://mikelove.github.io/fishpond/