Qui puoi trovare un codice R di esempio per calcolare la lunghezza del gene dato un file GTF (calcola anche il contenuto GC, che non ti serve). Questo utilizza uno dei numerosi modi per calcolare la lunghezza del gene, in questo caso la lunghezza del "modello di gene union". In questo metodo, gli esoni non duplicati per ciascun gene vengono semplicemente sommati ("non duplicati" in quanto nessuna base genomica viene conteggiata due volte). Questo è un modo molto semplice per ottenere una lunghezza del gene.

Esistono metodi alternativi di cui dovresti essere a conoscenza, tra i quali:

1. Lunghezza mediana della trascrizione: ovvero , le lunghezze esoniche in ciascuna trascrizione vengono sommate e viene utilizzata la mediana tra le trascrizioni. Questo è probabilmente un po 'più valido del codice a cui ho collegato.
2. Lunghezze geniche effettive per campione: il metodo ottimale, sebbene richieda l'uso di qualcosa come RSEM, che ti darà una lunghezza genica effettiva.
3. Usare la lunghezza dell '"isoforma maggiore" nel tessuto di interesse. Non è buono come il metodo 2, ma è più accurato di tutti gli altri.

Alla fine della giornata, stai solo elaborando un fattore di scala per ciascuno gene, quindi a meno che tu non intenda confrontare i valori tra i geni (questo è problematico per cominciare), è discutibile se usare alcuni dei metodi più corretti ma anche più impegnativi in ​​termini di tempo ti stia davvero ottenendo qualcosa.

Presumo che tu stia mappando contro il genoma piuttosto che il trascrittoma, poiché per il secondo la lunghezza sarebbe banale.

Supponendo il primo, penso che non dovrebbero essere incluse solo le sezioni di codifica ma anche l'UTR , dal momento che le letture possono essere mappate rispetto ad esse, che è ciò che alla fine ci interessa.

In generale, ho scoperto che i file di annotazione genica (ad esempio gff o gtf) possono essere incoerenti in termini di denominazione, quindi è buona pratica ispezionare e ricontrolla. Di seguito è riportato del codice R per importare l'annotazione e calcolare le lunghezze delle isoforme:

  library (rtracklayer) # Reading into a GRanges objectanno <- import.gff3 ("annotation.gff") # Filtering exons and UTRsexons <- anno [anno @ elementMetadata $ type% in% c ("exon", "five_prime_UTR", "three_prime_UTR"),]  

A seconda dell'annotazione a portata di mano, la più sensata è probabilmente meglio contare la lunghezza di ciascuna isoforma che è spesso contenuta nella colonna "Parent" del file di annotazione:

  # dividere le isoforme come preparazione per il passaggio successivo <- split (esoni , as.character (exons $ Parent)) # per ogni isoforma, calcola la somma di tutti gli esoni ridotti Gene_length <- sum (width (reduce (tmp)))  

Nota, riduci le unioni che si sovrappongono intervalli insieme, poiché gli UTR possono "contenere" bit di esoni che altrimenti verrebbero contati due volte.

Questo codice può ovviamente essere adattato principalmente cambiando "Parent", "exon" ecc.