Mercurial > repos > fubar > jbrowse2
view pafcount.py @ 22:2ddd41a0c2d5 draft
planemo upload for repository https://github.com/usegalaxy-eu/temporary-tools/tree/master/jbrowse2 commit be4f98b07024b59ff5e1ae0d8b467eecb15c7521-dirty
author | fubar |
---|---|
date | Thu, 01 Feb 2024 01:58:58 +0000 |
parents | |
children | 39b717d934a8 |
line wrap: on
line source
import sys """ Col Type Description 1 string Query sequence name 2 int Query sequence length 3 int Query start (0-based; BED-like; closed) 4 int Query end (0-based; BED-like; open) 5 char Relative strand: "+" or "-" 6 string Target sequence name 7 int Target sequence length 8 int Target start on original strand (0-based) 9 int Target end on original strand (0-based) 10 int Number of residue matches 11 int Alignment block length 12 int Mapping quality (0-255; 255 for missing) """ qcis = {} tcis = {} qtrans = {} ttrans = {} pafname = sys.argv[1] pf = open(pafname, 'r').readlines() for row in pf: qn,ql,qs,qe,qrs,tn,tl,ts,te,nm,abl,mq = row.strip().split("\t")[:12] if (qn == tn): # cis print('cis', qn,tn) tcis.setdefault(tn, 0) tcis[tn] = tcis[tn] + 1 qcis.setdefault(qn, 0) qcis[qn] = qcis[qn] + 1 else: # trans print('trans', qn,tn) k = '%s_%s' % (qn,tn) ttrans.setdefault(k, 0) ttrans[k] = ttrans[k]+ 1 qtrans.setdefault(k, 0) qtrans[k] = qtrans[k] + 1 #print('qcis', qcis,'\nqtrans', qtrans,'\ntcis', tcis,'\ntt', ttrans) #print('\nqtrans', qtrans,'\nttrans', ttrans) chroms = list(qtrans.keys()) print('chroms=', chroms) #print('chrom\tqcis\ttcis\tqtrans\tttrans') print('chrom\tqtrans\tttrans') for cr in chroms: #print('%s\t%d\t%d\t%d\t%d' % (cr, qcis[cr], tcis[cr], qtrans[cr], ttrans[cr])) print('%s\t%d\t%d' % (cr, qtrans[cr], ttrans[cr]))