Newer
Older

nina.marthe_ird.fr
committed
import argparse
import os
import subprocess
from getSegmentsCoordinates import seg_coord
def dir_path(path):
if os.path.isdir(path):
return path
else:
raise argparse.ArgumentTypeError(f"readable_dir:{path} is not a valid path")

nina.marthe_ird.fr
committed
def range_type(astr, min=0, max=100):
value = int(astr)
if min<= value <= max:
return value
else:
raise argparse.ArgumentTypeError('value not in range %s-%s'%(min,max))

nina.marthe_ird.fr
committed
def arg():
first_parser=argparse.ArgumentParser(add_help=False)
first_parser.add_argument('-pav','--pav_matrix',help="output a presence-absence variation matrix to recapitulate the annotation transfers; requires the -ann/--annotation option",action='store_true')

nina.marthe_ird.fr
committed
# annotation transfer options
first_parser.add_argument('-id', '--identity', metavar="[0-100]", type=range_type, choices=range(0,101), help="minimum identity percentage required for the feature to be transfered, between 0 and 100; default is 80; requires the -ann/--annotation option")
first_parser.add_argument('-cov', '--coverage', metavar="[0-100]", type=range_type, choices=range(0,101), help="minimum coverage percentage required for the feature to be transfered, between 0 and 100; default is 80; requires the -ann/--annotation option")
first_parser.add_argument('-diff', '--max_difference', metavar="N", type=int, default=0, help="maximum difference accepted for the feature to be transfered (transfer rejected if the feature size on the target genome is N times bigger or smaller than on the source genome; requires the -ann/--annotation option)")
first_args, _ = first_parser.parse_known_args()
parser = argparse.ArgumentParser(prog='GrAnnoT',description='Annotation transfer between genome and pangenome graph.',parents=[first_parser],usage='%(prog)s graph_file source_annotation_file source_genome_name [options]')

nina.marthe_ird.fr
committed
# input files
parser.add_argument('graph', metavar='graph.gfa', type=argparse.FileType('r'), help='pangenome graph file in GFA format containing the genome that the annotation refers to and the target genomes. the GFA file must have W-lines to describe the genomes walks in the graph, and not P-lines')

nina.marthe_ird.fr
committed
parser.add_argument('gff', metavar='annotation.gff', type=argparse.FileType('r') ,help='annotation file in GFF format containing the annotations to be transfered')
# optionnal input file
parser.add_argument('-coord', '--segment_coordinates_path', metavar="path/to/files", type=dir_path, default=".", help="path to the files given by the preprocessing of the graph; if not given the program will do the preprocessing; recommended if the program will run several times")
# genome names
parser.add_argument('source_genome', metavar="source_genome", type=str, help="name of the annotated genome's path in the GFA (field ? in the W line), that the annotation file refers to")
parser.add_argument('-t', '--target', metavar="target_genome", nargs='*', type=str, default=[], help="name of the target genomes' paths in the GFA (field ? in the W line) for the annotation; if the option is not used the programm will transfer the annotation on all the genomes")
# graph annotation options
parser.add_argument('-gff','--graph_gff',help="output the annotation of the graph in GFF format",action='store_true')
parser.add_argument('-gaf','--graph_gaf',help="output the annotation of the graph in GAF format",action='store_true')
# genome annotation options

nina.marthe_ird.fr
committed
parser.add_argument('-ann','--annotation',help="output the annotation transfer in GFF format",action='store_true',required=(first_args.pav_matrix or first_args.identity!=None or first_args.coverage!=None))

nina.marthe_ird.fr
committed
parser.add_argument('-var','--variation',help="output the detail of the variations for the transfered features",action='store_true')
parser.add_argument('-aln','--alignment',help="output the alignment of the transfered features",action='store_true')
# version
parser.add_argument('--version', action='version', version='%(prog)s 0.1')
return parser.parse_args()
def read_args(args):

nina.marthe_ird.fr
committed
if args.identity==None:
args.identity=80
if args.coverage==None:
args.coverage=80
path_list=args.target[:]

nina.marthe_ird.fr
committed
path_list.append(args.source_genome)
if args.segment_coordinates_path==".":
print("run getSegCoord.py")
seg_coord(args.graph.name,path_list)

nina.marthe_ird.fr
committed
elif args.segment_coordinates_path[-1]!="/":
args.segment_coordinates_path+="/"

nina.marthe_ird.fr
committed
print("run intersect") # intersect between gff and several bed
list_files=os.scandir(args.segment_coordinates_path)
command="echo -n "" > intersect"
subprocess.run(command,shell=True,timeout=None)
for file in list_files:
if args.source_genome in file.name:
command="bedtools intersect -wo -a "+file.path+" -b "+args.gff.name+">>intersect"
subprocess.run(command,shell=True,timeout=None)
#command="for f in "+args.segment_coordinates_path+"*"+args.source_genome+"*; do bedtools intersect -wo -a $f -b "+args.gff.name+" >> "+args.segment_coordinates_path+"intersect;done"
def get_genome_name(target_genomes,file_name):
for genome in target_genomes:
if genome in file_name:
return genome
return ""