%%bash
echo "TODAY'S DATE:"
date
echo "------------"
echo ""
#Display operating system info
lsb_release -a
echo ""
echo "------------"
echo "HOSTNAME: "; hostname
echo ""
echo "------------"
echo "Computer Specs:"
echo ""
lscpu
echo ""
echo "------------"
echo ""
echo "Memory Specs"
echo ""
free -mh
TODAY'S DATE:
Tue Dec 3 20:42:13 PST 2019
------------
Distributor ID: Ubuntu
Description: Ubuntu 16.04.6 LTS
Release: 16.04
Codename: xenial
------------
HOSTNAME:
swoose
------------
Computer Specs:
Architecture: x86_64
CPU op-mode(s): 32-bit, 64-bit
Byte Order: Little Endian
CPU(s): 24
On-line CPU(s) list: 0-23
Thread(s) per core: 2
Core(s) per socket: 6
Socket(s): 2
NUMA node(s): 1
Vendor ID: GenuineIntel
CPU family: 6
Model: 44
Model name: Intel(R) Xeon(R) CPU X5670 @ 2.93GHz
Stepping: 2
CPU MHz: 2926.194
BogoMIPS: 5852.50
Virtualization: VT-x
L1d cache: 32K
L1i cache: 32K
L2 cache: 256K
L3 cache: 12288K
NUMA node0 CPU(s): 0-23
Flags: fpu vme de pse tsc msr pae mce cx8 apic sep mtrr pge mca cmov pat pse36 clflush dts acpi mmx fxsr sse sse2 ss ht tm pbe syscall nx pdpe1gb rdtscp lm constant_tsc arch_perfmon pebs bts rep_good nopl xtopology nonstop_tsc aperfmperf pni pclmulqdq dtes64 monitor ds_cpl vmx smx est tm2 ssse3 cx16 xtpr pdcm pcid dca sse4_1 sse4_2 popcnt aes lahf_lm epb ssbd ibrs ibpb stibp kaiser tpr_shadow vnmi flexpriority ept vpid dtherm ida arat flush_l1d
------------
Memory Specs
total used free shared buff/cache available
Mem: 70G 4.4G 50G 559M 15G 65G
Swap: 4.7G 0B 4.7G
No LSB modules are available.
%env are best for bash
# Set workding directory
%env wd=/home/sam/analyses/20191203_pgen_v074.a4_genome_feature_counts
wd="/home/sam/analyses/20191203_pgen_v074.a4_genome_feature_counts"
# File download
%env downloads_list=Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3 Pgenerosa_v074.fa.fai
# Input/output files
%env old_fa_index=Pgenerosa_v074.fa.fai
%env org_merged_gff=Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3
%env new_merged_gff=Panopea-generosa-v1.0.a4.gff3
%env new_file=Panopea-generosa-v1.0.a4
%env new_repeats_file=Panopea-generosa-v1.0.a4.repeats
%env repeats_gff=Panopea-generosa-v1.0.a4.repeat_region.gff3
# Set lists of column header names
gff_header = ['seqid','source','type','start','end','score','strand','phase','attributes']
# Set genome size to 942Mbp
GENOME_SIZE = 942000000
env: wd=/home/sam/analyses/20191203_pgen_v074.a4_genome_feature_counts env: downloads_list=Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3 Pgenerosa_v074.fa.fai env: old_fa_index=Pgenerosa_v074.fa.fai env: org_merged_gff=Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3 env: new_merged_gff=Panopea-generosa-v1.0.a4.gff3 env: new_file=Panopea-generosa-v1.0.a4 env: new_repeats_file=Panopea-generosa-v1.0.a4.repeats env: repeats_gff=Panopea-generosa-v1.0.a4.repeat_region.gff3
import fnmatch
import os
import pandas
# Function to calculate percentage of genome comprised of a given feature
def ind_repeats_percent(feature_length_sum):
return round(float(feature_length_sum / GENOME_SIZE * 100), 2)
%%bash
mkdir --parents ${wd}
cd {wd}
/home/sam/analyses/20191203_pgen_v074.a4_genome_feature_counts
%%bash
mapfile -t downloads_array < <(echo ${downloads_list} | tr " " "\n")
for file in "${downloads_array[@]}"
do
rsync --archive --verbose \
owl:/volume1/web/halfshell/genomic-databank/${file} .
done
ls -lh
## If need to download via wget, uncomment and run lines below:
#for file in "${downloads_array[@]}"
#do
# wget --quiet --no-directories --no-check-certificate https://owl.fish.washington.edu/halfshell/genomic-databank/"${file}"
#done
receiving incremental file list Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3 sent 30 bytes received 580,544,162 bytes 20,369,971.65 bytes/sec total size is 580,473,160 speedup is 1.00 receiving incremental file list sent 11 bytes received 65 bytes 50.67 bytes/sec total size is 1,230 speedup is 16.18 total 554M -rw-rw-r-- 1 sam users 554M Dec 3 14:06 Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3 -rw-rw-rw- 1 sam users 1.3K Jun 26 08:54 Pgenerosa_v074.fa.fai
%%bash
# Array of old scaffold names
# Uses first column of FastA index file to get scaffold names
mapfile -t orig_scaffold_names < <(awk '{print $1}' "${old_fa_index}")
# Array of new scaffold names
mapfile -t new_scaffold_names < <(for number in {01..18}; do echo "Scaffold_${number}"; done)
# sed substitution
# creates sed script to find original scaffold names and replace them with new scafold names
# and passes to sed via stdin
for index in "${!orig_scaffold_names[@]}"
do
printf "s/%s/%s/\n" "${orig_scaffold_names[index]}" "${new_scaffold_names[index]}"
done | sed --file - "${org_merged_gff}" \
>> "${new_merged_gff}"
# Check that substituion worked
echo ""
echo "List of scaffold names in ${org_merged_gff}:"
echo "--------------------"
echo ""
awk '{print $1}' "${org_merged_gff}" | sort | uniq -c
echo ""
echo "--------------------"
echo ""
echo "List of scaffold names in ${new_merged_gff}:"
echo "--------------------"
echo ""
awk '{print $1}' "${new_merged_gff}" | sort | uniq -c
echo ""
List of scaffold names in Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3:
--------------------
1 ##Generated
1 ##gff-version
130305 PGA_scaffold10__49_contigs__length_53961475
118051 PGA_scaffold11__79_contigs__length_51449921
117911 PGA_scaffold12__71_contigs__length_50438331
109339 PGA_scaffold13__52_contigs__length_44396874
115751 PGA_scaffold14__91_contigs__length_45393038
111166 PGA_scaffold15__101_contigs__length_47938513
74528 PGA_scaffold16__33_contigs__length_31980953
84154 PGA_scaffold17__51_contigs__length_34923512
207214 PGA_scaffold1__77_contigs__length_89643857
115737 PGA_scaffold18__69_contigs__length_27737463
160589 PGA_scaffold2__36_contigs__length_69596280
136262 PGA_scaffold3__111_contigs__length_57743597
145359 PGA_scaffold4__129_contigs__length_65288255
150064 PGA_scaffold5__109_contigs__length_67248332
138842 PGA_scaffold6__104_contigs__length_61759565
101575 PGA_scaffold7__69_contigs__length_43120122
138958 PGA_scaffold8__63_contigs__length_61151155
84829 PGA_scaffold9__45_contigs__length_38581958
1 ##Project
--------------------
List of scaffold names in Panopea-generosa-v1.0.a4.gff3:
--------------------
1 ##Generated
1 ##gff-version
1 ##Project
207214 Scaffold_01
160589 Scaffold_02
136262 Scaffold_03
145359 Scaffold_04
150064 Scaffold_05
138842 Scaffold_06
101575 Scaffold_07
138958 Scaffold_08
84829 Scaffold_09
130305 Scaffold_10
118051 Scaffold_11
117911 Scaffold_12
109339 Scaffold_13
115751 Scaffold_14
111166 Scaffold_15
74528 Scaffold_16
84154 Scaffold_17
115737 Scaffold_18
%%bash
echo "Here are the features and their counts in ${new_merged_gff}:"
awk 'NR>3 {print $3}' "${new_merged_gff}" | sort | uniq -c
echo ""
# Create array for feature names
mapfile -t features_array < <(awk 'NR>3 {print $3}' "${new_merged_gff}" | sort | uniq)
# Save features in array
for feature in "${features_array[@]}"
do
output="${new_file}.${feature}.gff3"
head -n 3 "${new_merged_gff}" >> "${output}"
awk -v feature="$feature" '$3 == feature {print}' "${new_merged_gff}" \
>> "${output}"
echo ""
echo ""
echo "${output}"
echo "----------------------------------------------"
head -n 5 "${output}"
done
Here are the features and their counts in Panopea-generosa-v1.0.a4.gff3:
236960 CDS
236960 exon
34947 gene
38326 mRNA
1676544 repeat_region
8 rRNA
16889 tRNA
Panopea-generosa-v1.0.a4.CDS.gff3
----------------------------------------------
##gff-version 3
##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM
##Project Name : Pgenerosa_v074
Scaffold_01 GenSAS_5d9637f372b5d-publish CDS 2 125 . + 0 ID=PGEN_.00g000010.m01.CDS01;Name=PGEN_.00g000010.m01.CDS01;Parent=PGEN_.00g000010.m01;original_ID=cds.21510-PGEN_.00g234140.m01;Alias=cds.21510-PGEN_.00g234140.m01
Scaffold_01 GenSAS_5d9637f372b5d-publish CDS 1995 2095 . + 1 ID=PGEN_.00g000010.m01.CDS02;Name=PGEN_.00g000010.m01.CDS02;Parent=PGEN_.00g000010.m01;original_ID=cds.21510-PGEN_.00g234140.m01;Alias=cds.21510-PGEN_.00g234140.m01
Panopea-generosa-v1.0.a4.exon.gff3
----------------------------------------------
##gff-version 3
##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM
##Project Name : Pgenerosa_v074
Scaffold_01 GenSAS_5d9637f372b5d-publish exon 2 125 . + . ID=PGEN_.00g000010.m01.exon01;Name=PGEN_.00g000010.m01.exon01;Parent=PGEN_.00g000010.m01;original_ID=21510-PGEN_.00g234140.m01.exon1;Alias=21510-PGEN_.00g234140.m01.exon1
Scaffold_01 GenSAS_5d9637f372b5d-publish exon 1995 2095 . + . ID=PGEN_.00g000010.m01.exon02;Name=PGEN_.00g000010.m01.exon02;Parent=PGEN_.00g000010.m01;original_ID=21510-PGEN_.00g234140.m01.exon2;Alias=21510-PGEN_.00g234140.m01.exon2
Panopea-generosa-v1.0.a4.gene.gff3
----------------------------------------------
##gff-version 3
##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM
##Project Name : Pgenerosa_v074
Scaffold_01 GenSAS_5d9637f372b5d-publish gene 2 4719 . + . ID=PGEN_.00g000010;Name=PGEN_.00g000010;original_ID=21510-PGEN_.00g234140;Alias=21510-PGEN_.00g234140;original_name=21510-PGEN_.00g234140;Notes=sp|Q86IC9|CAMT1_DICDI [BLAST protein vs protein (blastp) 2.7.1],PF01596.12 [Pfam 1.6]
Scaffold_01 GenSAS_5d9637f372b5d-publish gene 19808 36739 . - . ID=PGEN_.00g000020;Name=PGEN_.00g000020;original_ID=21510-PGEN_.00g234150;Alias=21510-PGEN_.00g234150;original_name=21510-PGEN_.00g234150;Notes=sp|P04177|TY3H_RAT [BLAST protein vs protein (blastp) 2.7.1],sp|P04177|TY3H_RAT [DIAMOND Functional 0.9.22],IPR036951 [InterProScan 5.29-68.0],PF00351.16 [Pfam 1.6]
Panopea-generosa-v1.0.a4.mRNA.gff3
----------------------------------------------
##gff-version 3
##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM
##Project Name : Pgenerosa_v074
Scaffold_01 GenSAS_5d9637f372b5d-publish mRNA 2 4719 . + . ID=PGEN_.00g000010.m01;Name=PGEN_.00g000010.m01;Parent=PGEN_.00g000010;original_ID=21510-PGEN_.00g234140.m01;Alias=21510-PGEN_.00g234140.m01;original_name=21510-PGEN_.00g234140
Scaffold_01 GenSAS_5d9637f372b5d-publish mRNA 19808 36739 . - . ID=PGEN_.00g000020.m01;Name=PGEN_.00g000020.m01;Parent=PGEN_.00g000020;original_ID=21510-PGEN_.00g234150.m01;Alias=21510-PGEN_.00g234150.m01;original_name=21510-PGEN_.00g234150
Panopea-generosa-v1.0.a4.repeat_region.gff3
----------------------------------------------
##gff-version 3
##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM
##Project Name : Pgenerosa_v074
Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 265 338 334 + . ID=19535.GS22252505.PGEN_.repeat00000002;Name=19535.GS22252505.PGEN_.repeat00000002;repeat_match=rnd-5_family-367;repeat_class=Unknown;
Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 855 1050 451 + . ID=19535.GS22252505.PGEN_.repeat00000003;Name=19535.GS22252505.PGEN_.repeat00000003;repeat_match=rnd-5_family-1818;repeat_class=Unknown;
Panopea-generosa-v1.0.a4.rRNA.gff3
----------------------------------------------
##gff-version 3
##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM
##Project Name : Pgenerosa_v074
Scaffold_01 RNAmmer-1.2 rRNA 46743424 46743531 48.2 - . ID=21514.GS22252505.PGEN_.rRNA00000007;Name=21514.GS22252505.PGEN_.rRNA00000007;molecule_type=8s_rRNA;
Scaffold_02 RNAmmer-1.2 rRNA 69261970 69262107 51.4 + . ID=21514.GS22252506.PGEN_.rRNA00000003;Name=21514.GS22252506.PGEN_.rRNA00000003;molecule_type=8s_rRNA;
Panopea-generosa-v1.0.a4.tRNA.gff3
----------------------------------------------
##gff-version 3
##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM
##Project Name : Pgenerosa_v074
Scaffold_01 GenSAS_5d82b316cd298-trnascan tRNA 8317 8393 34.9 . . ID=21513.GS22252505.PGEN_.tRNA00000001;Name=Ala;anti_codon=CGC;gene_no=1;
Scaffold_01 GenSAS_5d82b316cd298-trnascan tRNA 87519 87602 34.9 . . ID=21513.GS22252505.PGEN_.tRNA00000002;Name=Thr;anti_codon=AGT;gene_no=2;
for file in os.listdir('.'):
if fnmatch.fnmatch(file, 'Panopea-generosa-v1.0.a4*.gff3'):
print('\n' * 2)
print(file)
print("-------------------------")
# Import GFF.
# Skip first 3 rows (gff header lines) and indicate file is tab-separated
gff=pandas.read_csv(file, header=None, skiprows=3, sep="\t")
# Rename columns
gff.columns = gff_header
# Subtract start value from end value.
# Have to add 1 so that sequence length can't equal zero (i.e. adjust for 1-based counting system)
gff['seqlength'] = gff.apply(lambda position: position['end'] - position['start'] + 1, axis=1)
# Apply functions in list to seqlength column
gff_stats = gff['seqlength'].agg(['mean', 'min', 'median', 'max'])
print (gff_stats)
Panopea-generosa-v1.0.a4.gene.gff3 ------------------------- mean 10811.04461 min 166.00000 median 4464.00000 max 283066.00000 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.mRNA.gff3 ------------------------- mean 12903.649559 min 166.000000 median 5453.000000 max 283066.000000 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.tRNA.gff3 ------------------------- mean 74.807745 min 53.000000 median 75.000000 max 314.000000 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.rRNA.gff3 ------------------------- mean 117.125 min 108.000 median 115.000 max 138.000 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.repeat_region.gff3 -------------------------
/home/sam/programs/minicocnda3/lib/python3.6/site-packages/IPython/core/interactiveshell.py:3049: DtypeWarning: Columns (5) have mixed types. Specify dtype option on import or set low_memory=False. interactivity=interactivity, compiler=compiler, result=result)
mean 212.244974 min 6.000000 median 149.000000 max 10981.000000 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.exon.gff3 ------------------------- mean 201.476988 min 3.000000 median 133.000000 max 13221.000000 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.gff3 ------------------------- mean 591.326038 min 3.000000 median 148.000000 max 283066.000000 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.CDS.gff3 ------------------------- mean 201.476988 min 3.000000 median 133.000000 max 13221.000000 Name: seqlength, dtype: float64
%%bash
# Initialize array
features_array=()
# Identify unique features in GFF
## Store as an array
## Skip first three header lines, and then cut on two delimiters that are present
echo "Unique repeats features in ${repeats_gff}:"
while IFS='' read -r line
do
features_array+=("$line")
done < <(awk -F"class=" 'NR >3 {print $2}' "${repeats_gff}" \
| sort -u \
| cut -d '%' -f 1 \
| cut -d ';' -f 1 \
| uniq)
# Check array contents
for feature in "${features_array[@]}"
do
echo "${feature}"
done
echo ""
echo "----------------------------------------------"
echo ""
# Loop through array and create new GFFs from each feature
for feature in "${features_array[@]}"
do
echo "Parsing ${feature} from ${repeats_gff}..."
echo "Writing GFF3 header to ${new_repeats_file}.${feature}.gff3"
# Write header to new file
head -n 3 "${repeats_gff}" > "${new_repeats_file}"."${feature}".gff3
echo "Parsing matching feature lines for ${feature} feature..."
grep "repeat_class=${feature}" "${repeats_gff}" >> "${new_repeats_file}"."${feature}".gff3
echo "Done with parsing ${feature} feature."
# Count number of lines, excluding three line header (oddly, need to run tail --lines n+1, as +1 prints entire file)
feature_count=$(tail --lines +4 "${new_repeats_file}"."${feature}".gff3 | wc -l)
echo "Identified ${feature_count} ${feature} features."
echo "Output file is: ${new_repeats_file}.${feature}.gff3"
head -n 5 "${new_repeats_file}.${feature}.gff3"
echo ""
echo "----------------------------------------------"
done
echo "----------------------------------------------"
echo ""
ls -lh
Unique repeats features in Panopea-generosa-v1.0.a4.repeat_region.gff3: DNA LINE LTR RC Simple_repeat SINE Unknown ---------------------------------------------- Parsing DNA from Panopea-generosa-v1.0.a4.repeat_region.gff3... Writing GFF3 header to Panopea-generosa-v1.0.a4.repeats.DNA.gff3 Parsing matching feature lines for DNA feature... Done with parsing DNA feature. Identified 23195 DNA features. Output file is: Panopea-generosa-v1.0.a4.repeats.DNA.gff3 ##gff-version 3 ##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM ##Project Name : Pgenerosa_v074 Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 11381 11673 2176 + . ID=19535.GS22252505.PGEN_.repeat00000029;Name=19535.GS22252505.PGEN_.repeat00000029;repeat_match=rnd-1_family-791;repeat_class=DNA%2FTcMar-Mariner; Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 29272 29508 341 - . ID=19535.GS22252505.PGEN_.repeat00000058;Name=19535.GS22252505.PGEN_.repeat00000058;repeat_match=rnd-1_family-443;repeat_class=DNA%2FTcMar-Mariner; ---------------------------------------------- Parsing LINE from Panopea-generosa-v1.0.a4.repeat_region.gff3... Writing GFF3 header to Panopea-generosa-v1.0.a4.repeats.LINE.gff3 Parsing matching feature lines for LINE feature... Done with parsing LINE feature. Identified 75939 LINE features. Output file is: Panopea-generosa-v1.0.a4.repeats.LINE.gff3 ##gff-version 3 ##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM ##Project Name : Pgenerosa_v074 Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 5885 6044 316 + . ID=19535.GS22252505.PGEN_.repeat00000014;Name=19535.GS22252505.PGEN_.repeat00000014;repeat_match=rnd-5_family-963;repeat_class=LINE%2FRTE-X; Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 45583 45756 460 - . ID=19535.GS22252505.PGEN_.repeat00000074;Name=19535.GS22252505.PGEN_.repeat00000074;repeat_match=rnd-5_family-33;repeat_class=LINE%2FL2; ---------------------------------------------- Parsing LTR from Panopea-generosa-v1.0.a4.repeat_region.gff3... Writing GFF3 header to Panopea-generosa-v1.0.a4.repeats.LTR.gff3 Parsing matching feature lines for LTR feature... Done with parsing LTR feature. Identified 3255 LTR features. Output file is: Panopea-generosa-v1.0.a4.repeats.LTR.gff3 ##gff-version 3 ##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM ##Project Name : Pgenerosa_v074 Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 1275 1466 920 - . ID=19535.GS22252505.PGEN_.repeat00000005;Name=19535.GS22252505.PGEN_.repeat00000005;repeat_match=rnd-5_family-1106;repeat_class=LTR%2FGypsy; Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 670854 671005 266 + . ID=19535.GS22252505.PGEN_.repeat00001290;Name=19535.GS22252505.PGEN_.repeat00001290;repeat_match=rnd-5_family-1106;repeat_class=LTR%2FGypsy; ---------------------------------------------- Parsing RC from Panopea-generosa-v1.0.a4.repeat_region.gff3... Writing GFF3 header to Panopea-generosa-v1.0.a4.repeats.RC.gff3 Parsing matching feature lines for RC feature... Done with parsing RC feature. Identified 603 RC features. Output file is: Panopea-generosa-v1.0.a4.repeats.RC.gff3 ##gff-version 3 ##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM ##Project Name : Pgenerosa_v074 Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 2890875 2891508 5041 - . ID=19535.GS22252505.PGEN_.repeat00005106;Name=19535.GS22252505.PGEN_.repeat00005106;repeat_match=rnd-1_family-384;repeat_class=RC%2FHelitron; Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 3120337 3120414 311 + . ID=19535.GS22252505.PGEN_.repeat00005530;Name=19535.GS22252505.PGEN_.repeat00005530;repeat_match=rnd-1_family-384;repeat_class=RC%2FHelitron; ---------------------------------------------- Parsing Simple_repeat from Panopea-generosa-v1.0.a4.repeat_region.gff3... Writing GFF3 header to Panopea-generosa-v1.0.a4.repeats.Simple_repeat.gff3 Parsing matching feature lines for Simple_repeat feature... Done with parsing Simple_repeat feature. Identified 19865 Simple_repeat features. Output file is: Panopea-generosa-v1.0.a4.repeats.Simple_repeat.gff3 ##gff-version 3 ##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM ##Project Name : Pgenerosa_v074 Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 8890 9134 853 + . ID=19535.GS22252505.PGEN_.repeat00000022;Name=19535.GS22252505.PGEN_.repeat00000022;repeat_match=rnd-1_family-158;repeat_class=Simple_repeat; Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 99891 100132 1511 + . ID=19535.GS22252505.PGEN_.repeat00000197;Name=19535.GS22252505.PGEN_.repeat00000197;repeat_match=rnd-6_family-10;repeat_class=Simple_repeat; ---------------------------------------------- Parsing SINE from Panopea-generosa-v1.0.a4.repeat_region.gff3... Writing GFF3 header to Panopea-generosa-v1.0.a4.repeats.SINE.gff3 Parsing matching feature lines for SINE feature... Done with parsing SINE feature. Identified 43129 SINE features. Output file is: Panopea-generosa-v1.0.a4.repeats.SINE.gff3 ##gff-version 3 ##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM ##Project Name : Pgenerosa_v074 Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 2928 3030 668 - . ID=19535.GS22252505.PGEN_.repeat00000010;Name=19535.GS22252505.PGEN_.repeat00000010;repeat_match=rnd-3_family-517;repeat_class=SINE%2FtRNA-Core-L2; Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 3050 3146 311 - . ID=19535.GS22252505.PGEN_.repeat00000011;Name=19535.GS22252505.PGEN_.repeat00000011;repeat_match=rnd-1_family-48;repeat_class=SINE%2FtRNA-Core-L2; ---------------------------------------------- Parsing Unknown from Panopea-generosa-v1.0.a4.repeat_region.gff3... Writing GFF3 header to Panopea-generosa-v1.0.a4.repeats.Unknown.gff3 Parsing matching feature lines for Unknown feature... Done with parsing Unknown feature. Identified 1510558 Unknown features. Output file is: Panopea-generosa-v1.0.a4.repeats.Unknown.gff3 ##gff-version 3 ##Generated using GenSAS, Tuesday 26th of November 2019 07:12:25 PM ##Project Name : Pgenerosa_v074 Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 265 338 334 + . ID=19535.GS22252505.PGEN_.repeat00000002;Name=19535.GS22252505.PGEN_.repeat00000002;repeat_match=rnd-5_family-367;repeat_class=Unknown; Scaffold_01 GenSAS_5d25089d78791-repeatmodeler repeat_region 855 1050 451 + . ID=19535.GS22252505.PGEN_.repeat00000003;Name=19535.GS22252505.PGEN_.repeat00000003;repeat_match=rnd-5_family-1818;repeat_class=Unknown; ---------------------------------------------- ---------------------------------------------- total 1.9G -rw-rw-r-- 1 sam sam 53M Dec 4 13:14 Panopea-generosa-v1.0.a4.CDS.gff3 -rw-rw-r-- 1 sam sam 55M Dec 4 13:14 Panopea-generosa-v1.0.a4.exon.gff3 -rw-rw-r-- 1 sam sam 9.5M Dec 4 13:14 Panopea-generosa-v1.0.a4.gene.gff3 -rw-rw-r-- 1 sam sam 486M Dec 4 13:13 Panopea-generosa-v1.0.a4.gff3 -rw-rw-r-- 1 sam sam 9.1M Dec 4 13:14 Panopea-generosa-v1.0.a4.mRNA.gff3 -rw-rw-r-- 1 sam sam 358M Dec 4 13:14 Panopea-generosa-v1.0.a4.repeat_region.gff3 -rw-rw-r-- 1 sam sam 5.2M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.DNA.gff3 -rw-rw-r-- 1 sam sam 17M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.LINE.gff3 -rw-rw-r-- 1 sam sam 726K Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.LTR.gff3 -rw-rw-r-- 1 sam sam 136K Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.RC.gff3 -rw-rw-r-- 1 sam sam 4.4M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.Simple_repeat.gff3 -rw-rw-r-- 1 sam sam 9.7M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.SINE.gff3 -rw-rw-r-- 1 sam sam 322M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.Unknown.gff3 -rw-rw-r-- 1 sam sam 1.4K Dec 4 13:15 Panopea-generosa-v1.0.a4.rRNA.gff3 -rw-rw-r-- 1 sam sam 2.5M Dec 4 13:15 Panopea-generosa-v1.0.a4.tRNA.gff3 -rw-rw-r-- 1 sam users 554M Dec 3 14:06 Panopea-generosa-vv0.74.a4-merged-2019-11-26-15-30-34.gff3 -rw-rw-rw- 1 sam users 1.3K Jun 26 08:54 Pgenerosa_v074.fa.fai
total_repeats_percent = 0
for file in os.listdir('.'):
if fnmatch.fnmatch(file, 'Panopea-generosa-v1.0.a4.repeats*.gff3'):
print('\n' * 2)
print(file)
print("-------------------------")
# Import GFF.
# Skip first five rows and file is tab-separated
gff=pandas.read_csv(file, header=None, skiprows=5, sep="\t")
# Rename columns
gff.columns = gff_header
# Subtract start value from end value.
# Have to add 1 so that sequence length can't equal zero
gff['seqlength'] = gff.apply(lambda position: position['end'] - position['start'] + 1, axis=1)
gff_sum = gff['seqlength'].sum()
total_repeats_percent += ind_repeats_percent(gff_sum)
print ("percent" , ind_repeats_percent(gff_sum))
# Apply functions in list to seqlength column
gff_stats = gff['seqlength'].agg(['sum', 'mean', 'min', 'median', 'max'])
print (gff_stats.round(2))
print('\n' * 2)
print("-------------------------")
print ("Repeats composition of genome (percent):" , total_repeats_percent)
Panopea-generosa-v1.0.a4.repeats.LTR.gff3 ------------------------- percent 0.25 sum 2315583.00 mean 711.83 min 11.00 median 323.00 max 6541.00 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.repeats.RC.gff3 ------------------------- percent 0.03 sum 258182.00 mean 429.59 min 13.00 median 464.00 max 674.00 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.repeats.Simple_repeat.gff3 ------------------------- percent 0.55 sum 5138701.00 mean 258.71 min 6.00 median 124.00 max 5981.00 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.repeats.DNA.gff3 ------------------------- percent 1.01 sum 9497156.00 mean 409.48 min 11.00 median 248.00 max 7012.00 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.repeats.SINE.gff3 ------------------------- percent 0.72 sum 6737909.00 mean 156.23 min 11.00 median 165.00 max 934.00 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.repeats.LINE.gff3 ------------------------- percent 3.19 sum 30035624.00 mean 395.53 min 11.00 median 226.00 max 6604.00 Name: seqlength, dtype: float64 Panopea-generosa-v1.0.a4.repeats.Unknown.gff3 ------------------------- percent 32.04 sum 3.018520e+08 mean 1.998300e+02 min 1.100000e+01 median 1.450000e+02 max 1.098100e+04 Name: seqlength, dtype: float64 ------------------------- Repeats composition of genome (percent): 37.79
%%bash
%%bash
mapfile -t downloads_array < <(echo ${downloads_list} | tr " " "\n")
for file in "${downloads_array[@]}"
do
rm "${file}"
done
ls -lh
total 1.3G -rw-rw-r-- 1 sam sam 53M Dec 4 13:14 Panopea-generosa-v1.0.a4.CDS.gff3 -rw-rw-r-- 1 sam sam 55M Dec 4 13:14 Panopea-generosa-v1.0.a4.exon.gff3 -rw-rw-r-- 1 sam sam 9.5M Dec 4 13:14 Panopea-generosa-v1.0.a4.gene.gff3 -rw-rw-r-- 1 sam sam 486M Dec 4 13:13 Panopea-generosa-v1.0.a4.gff3 -rw-rw-r-- 1 sam sam 9.1M Dec 4 13:14 Panopea-generosa-v1.0.a4.mRNA.gff3 -rw-rw-r-- 1 sam sam 358M Dec 4 13:14 Panopea-generosa-v1.0.a4.repeat_region.gff3 -rw-rw-r-- 1 sam sam 5.2M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.DNA.gff3 -rw-rw-r-- 1 sam sam 17M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.LINE.gff3 -rw-rw-r-- 1 sam sam 726K Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.LTR.gff3 -rw-rw-r-- 1 sam sam 136K Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.RC.gff3 -rw-rw-r-- 1 sam sam 4.4M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.Simple_repeat.gff3 -rw-rw-r-- 1 sam sam 9.7M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.SINE.gff3 -rw-rw-r-- 1 sam sam 322M Dec 4 13:19 Panopea-generosa-v1.0.a4.repeats.Unknown.gff3 -rw-rw-r-- 1 sam sam 1.4K Dec 4 13:15 Panopea-generosa-v1.0.a4.rRNA.gff3 -rw-rw-r-- 1 sam sam 2.5M Dec 4 13:15 Panopea-generosa-v1.0.a4.tRNA.gff3
bash: line 1: fg: no job control
%%bash
cd ..
rsync --progress --archive --relative \
./20191203_pgen_v074.a4_genome_feature_counts \
gannet:/volume2/web/Atumefaciens
sending incremental file list
./
20191203_pgen_v074.a4_genome_feature_counts/
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.CDS.gff3
55,346,545 100% 110.59MB/s 0:00:00 (xfr#1, to-chk=14/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.exon.gff3
57,148,143 100% 55.44MB/s 0:00:00 (xfr#2, to-chk=13/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.gene.gff3
9,888,462 100% 8.81MB/s 0:00:01 (xfr#3, to-chk=12/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.gff3
509,354,871 100% 105.46MB/s 0:00:04 (xfr#4, to-chk=11/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.mRNA.gff3
9,450,148 100% 12.99MB/s 0:00:00 (xfr#5, to-chk=10/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.rRNA.gff3
1,388 100% 1.95kB/s 0:00:00 (xfr#6, to-chk=9/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeat_region.gff3
374,972,067 100% 87.71MB/s 0:00:04 (xfr#7, to-chk=8/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeats.DNA.gff3
5,376,409 100% 39.75MB/s 0:00:00 (xfr#8, to-chk=7/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeats.LINE.gff3
17,343,693 100% 41.56MB/s 0:00:00 (xfr#9, to-chk=6/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeats.LTR.gff3
743,047 100% 1.75MB/s 0:00:00 (xfr#10, to-chk=5/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeats.RC.gff3
138,457 100% 333.03kB/s 0:00:00 (xfr#11, to-chk=4/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeats.SINE.gff3
10,113,860 100% 19.45MB/s 0:00:00 (xfr#12, to-chk=3/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeats.Simple_repeat.gff3
4,538,129 100% 8.07MB/s 0:00:00 (xfr#13, to-chk=2/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.repeats.Unknown.gff3
336,719,168 100% 75.65MB/s 0:00:04 (xfr#14, to-chk=1/17)
20191203_pgen_v074.a4_genome_feature_counts/Panopea-generosa-v1.0.a4.tRNA.gff3
2,548,814 100% 10.48MB/s 0:00:00 (xfr#15, to-chk=0/17)