{ Les snippets de bioinfo-fr.net }


Parseur de fichier gff3 basique [Python]

08.07.2018     propan2one      gff3 parser python3 

  La fonction gff3parser() permet de parser un fichier au format .gff3 et de retourner un dictionnaire au format :

{ numéro de l'ORF : [start, stop, strand ], ... }

Les print en commentaire permettent d'afficher le contenue du gff3 sur le terminal

contentGFF = ["seqid", "source", "type", "start", "end", "score", "strand", "phase", "attributes"]
def gff3parser (filename, dictOrf):
	 """ parser de gff3 """
	with open(filename) as f:
		for line in f:
			if line.startswith("#"):
		                continue
		content = line.strip().split("\t")
		if len(content) == len(contentGFF):
			value = [content[3], content[4], content[6]]
			#print(value)
			#print(content[8])
			numORF = content[8].strip().split(";")
			#print(numORF)
			for product in numORF:
				if re.search(r"^product=ORF\d{0,3}", product):
					keys = product[11:]
					dictOrf[keys] = value

		return {}

# Usage :
orf = {}
gff3parser(gffname,orf) #avec gffname le nom du fichier gff3 à parser
1/5 - [1 vote]




Retire le numéro de version du code Ensembl des gènes / transcrits [R]

05.07.2018     gdevailly      GENCODE Ensembl 

  Les codes Ensembl des gènes / transcrits viennent parfois avec un chiffre final traçant le nombre de modifications de ces annotations:

Par exemple, l'annotation de ENSG00000134046.11 a été modifié 11 fois.

Cette petit expression régulière permet de s'en débarrasser, pour faciliter les croisements de listes de gènes.
x <- c("ENSG00000134046.11", "ENSG00000229807.10")

sub("\\.[0-9]*$", "", x)
# [1] "ENSG00000134046" "ENSG00000229807"
0/5 - [0 vote]




Parseur de csv basique [Python]

04.07.2018     qcavaille      python3 parser csv 

  Parseur de csv basique :
- Délimiteur réglable
- Pas de limite de colonnes, mais nécessite des lignes "complètes" ET une ligne d'en-têtes de colonnes
- Retourne un générateur, théoriquement plus léger en mémoire. Mais itérable une seule fois !
import csv as L
def readCSV(file):
    with open(file) as csvfile:
        readCSV = L.reader(csvfile, delimiter=";")
        firstRow=True
        columnName=[]
        for row in readCSV:
            if firstRow:
                firstRow=False
                for col in row:
                    columnName.append(col)
            else:
                aRow={}
                i=0
                for col in columnName:
                    aRow[col]=row[i]
                    i+=1
                yield aRow
2/5 - [3 votes]




reverse complement a sequence (probably the best solution or near of them) [C++]

04.07.2018     Natir      rev_comp singleton 

  A powerful reverse complement algorithm, encapsulating it in a singleton to avoid unnecessary multiple allocation.

#include <string>

class rev_comp
{
public:
    static std::string run(std::string seq)
    {
	rev_comp::build_instance();
	
	auto first = seq.begin(), last = seq.end();
        
	while(true)
        {
		if(first == last || first == --last)
        	{
        	    if(seq.length() % 2)
        		    *first = rev_comp::complement[*first];
        	    return seq;
        	}
        	else
        	{
        	    *first = rev_comp::complement[*first];
        	    *last = rev_comp::complement[*last];
        	    std::iter_swap(first, last);
        	    ++first;
            }
        }
    }

protected:
    static rev_comp* _instance;
    static void build_instance()
    {
        if(_instance == nullptr)
            _instance = new rev_comp();
    }
    
    rev_comp()
    {
        this->complement['A'] = 'T';
        this->complement['T'] = 'A';
        this->complement['C'] = 'G';
        this->complement['G'] = 'C';
        this->complement['a'] = 't';
        this->complement['t'] = 'a';
        this->complement['c'] = 'g';
        this->complement['g'] = 'c';
    }
        
    ~rev_comp(){}
    
    
private:
    static char complement['t' + 1];
};

char rev_comp::complement['t' + 1];
rev_comp* rev_comp::_instance = nullptr;
5/5 - [2 votes]




Diminution de la taille d'une matrice [R]

04.07.2018     gdevailly      matrice interpolation heatmap 

  Ploter une grosse matrice sous forme d'heatmap peut prendre du temps.
Il peux être judicieux de réduire la taille de la matrice avant de la ploter, en moyennant des groupes de cellules adjacentes (ou en prenant le maximum, lorsque la matrice est 'sparse').
C'est notamment utile lorsqu'il y a beaucoup plus de cellules que de pixels sur les écran actuels (par exemple, ploter 80 000 régions génomiques sur un écran HD de 1080 pixels de hauts).

Utilisation:

> bigmat <- matrix(1:64, ncol = 8)
> bigmat
[,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8]
[1,] 1 9 17 25 33 41 49 57
[2,] 2 10 18 26 34 42 50 58
[3,] 3 11 19 27 35 43 51 59
[4,] 4 12 20 28 36 44 52 60
[5,] 5 13 21 29 37 45 53 61
[6,] 6 14 22 30 38 46 54 62
[7,] 7 15 23 31 39 47 55 63
[8,] 8 16 24 32 40 48 56 64
> # default usage
> redim_matrix(bigmat, target_height = 4, target_width = 3)
[,1] [,2] [,3]
[1,] 10.0 30.0 50.0
[2,] 11.5 31.5 51.5
[3,] 13.0 33.0 53.0
[4,] 15.0 35.0 55.0
> # changing aggregating function
> redim_matrix(bigmat, target_height = 4, target_width = 3, summary_func = function(x) max(x, na.rm = TRUE))
[,1] [,2] [,3]
[1,] 19 43 59
[2,] 20 44 60
[3,] 22 46 62
[4,] 24 48 64
> # multicore
> redim_matrix(bigmat, target_height = 4, target_width = 3, n_core = 2)
[,1] [,2] [,3]
[1,] 10.0 30.0 50.0
[2,] 11.5 31.5 51.5
[3,] 13.0 33.0 53.0
[4,] 15.0 35.0 55.0
# reduce matrix size, using a summarizing function (default, mean)
redim_matrix <- function(
    mat,
    target_height = 100,
    target_width = 100,
    summary_func = function(x) mean(x, na.rm = TRUE),
    n_core = 1
    ) {

    if(target_height > nrow(mat) | target_width > ncol(mat)) {
        stop("Input matrix must be bigger than target width and height.")
    }

    seq_height <- round(seq(1, nrow(mat), length.out = target_height + 1))
    seq_width  <- round(seq(1, ncol(mat), length.out = target_width  + 1))

    # complicate way to write a double for loop
    do.call(rbind, parallel::mclapply(seq_len(target_height), function(i) { # i is row
        vapply(seq_len(target_width), function(j) { # j is column
            summary_func(
                mat[
                    seq(seq_height[i], seq_height[i + 1]),
                    seq(seq_width[j] , seq_width[j + 1] )
                    ]
            )
        }, 0.0)
    }, mc.cores = n_core))
}
5/5 - [1 vote]




Réarrangement de l'ordre des branches d'un dendrogramme [R]

04.07.2018     gdevailly      hclust dendrogram optimal leaf ordering OLO seriation clustering heatmap 

  Lors d'un clustering hierarchique, l'agencement des branches à chaque nœud est laissé libre (telle branche peut être agencé à gauche ou à droite).
Il peut en résulter des visualisations de heatmaps clusterisée parfois pas très belle.
Heureusement, il est possible de réarranger les branches d'un clustering hiérarchique pour aboutir à de jolies images, tel qu'illustré dans ce billet de blog:
http://nicolas.kruchten.com/content/2018/02/seriation/

Cette petite fonction remplacera donc avantageusement la fonction hclust() dans vos heatmaps clusterisées.
# perform hclust, then aply optimal leaf ordering of the dendrogram branches
require(seriation) # you should install.packages("seriation") if error message

hclust_olo <- function(mdist, ...) {
    myClust <- hclust(mdist, ...)
    myOlo <- seriation::seriate(mdist, method = "OLO")
    seriation::permute(myClust, myOlo)
}
5/5 - [1 vote]





<<< Snippets récents
Powered by canSnippet CE - by ademcan