php 统计fasta 序列长度和GC含量

最近php7的消息铺天盖地, 忍不住想尝试下。星期天看了下语法, 写个小脚本练下手;

这个脚本读取fasta 文件, 输出序列的长度和GC含量;

<?php
$fasta = "test.fasta";
$meta  = array();
$meta  = parse_fasta($fasta);
write_res($meta);

function parse_fasta($fasta) {
    $meta = array();
    $file_handle = fopen($fasta, 'r');
    $id = '';
    $seq = '';
    while (!feof($file_handle)) {
        $line = fgets($file_handle);
        $line = preg_replace("/\s+$/", "", $line);
        if (preg_match("/^>/", $line)) {
            if ($id) {
                $meta[$id] = $seq;
                $id = $line;
                $seq = '';
            } else {
                $id = $line;
        }     
        } else {
            $seq .= $line;
        }
    }
    $meta[$id] = $seq;
    fclose($file_handle);
    return $meta;
}


function write_res($meta) {
    foreach ($meta as $key => $value) {
        $len = cal_length($value);
        $gc  = cal_gc($value);
        echo "$key\t$value\t$len\t$gc\n";
    }
}


function cal_length($seq) {
    return strlen($seq);

}


function cal_gc($seq) {
    $gc = array();
    preg_match_all("/G|C|g|c/", $seq, $gc);
    return count($gc[0]) / strlen($seq);
}


?>

 

posted on 2015-12-08 12:24  庐州月光  阅读(751)  评论(0编辑  收藏  举报