We have the fasta files of some viral genomes. We want to aligend all of them and extract only one gene of them.
cat *.fasta > all.fasta
cat gene.fasta all.fasta > all_with_gene.fasta
The I ran mafft on them
downlaod and unzip htslib | |
/Users/sina/Downloads/htslib-1.10.2.tar.bz2 | |
cd /Users/sina/Downloads/htslib-1.10.2 | |
mkdir /installation_folder | |
./configure --prefix=/Users/sina/Downloads/htslib-1.10.2/installation_folder | |
make |
We have the fasta files of some viral genomes. We want to aligend all of them and extract only one gene of them.
cat *.fasta > all.fasta
cat gene.fasta all.fasta > all_with_gene.fasta
The I ran mafft on them
#!/usr/bin/env python3 | |
"""Coverage mean and standard deviation of autosomes | |
Estimate the mean and standard deviation from a mosdepth coverage BED for | |
positions with coverage in the range (0, 2 * non-zero mode). This estimate | |
behaves well for PacBio HiFi WGS of human germline aligned to either hs37d5 and | |
GRCh38, and may be useful for other situations as well. | |
$ bash mosdepth --threads 3 --no-per-base --by 500 -m "${BAM%.*}.median" "${BAM}" | |
$ tabix ${BAM%.*}.median.regions.bed.gz {1..22} | python depth_mean_stddev.py |
#!/usr/bin/python3 | |
import numpy as np | |
from sys import argv | |
file_fq_input_addrss = argv[1] | |
file_fq_output_addrss = argv[2] | |
file_fq_input= open(file_fq_input_addrss,'r'); |
A simpe code for extracting the regions with high coverage (>depth_treshold
) from a bed file.
The input of the code is a bed file containing rows of regions with its depth as following format
21 9424000 9424500 69.00
21 9424500 9425000 69.00
21 9425000 9425500 67.00
21 9425500 9426000 66.00
21 9426000 9426500 65.00
21 9426500 9427000 66.00