Commit 171a1aed authored by fabrice's avatar fabrice
Browse files

Bug Fixe: N_factor convert to int during mapfilter step. Float values wasn't...

Bug Fixe: N_factor convert to int during mapfilter step. Float values wasn't taken into account before ...
parent 02494b4e
......@@ -195,7 +195,7 @@ class AriaEcCommand:
default=False,
help="Use secondary structure index")
parser.add_argument("--prefix", dest="prefix", default="",
help="File name")
help="Contact map name", nargs="+")
return parser
def create_settings(self):
......
No preview for this file type
......@@ -131,7 +131,6 @@ pickle_output: no
[contactmap]
; -------------------------- Contactmap parameters --------------------------- #
; Plot settings
n_factor: 1.5
save_fig: True
size_fig: 10
plot_ext: pdf
......
......@@ -37,13 +37,13 @@ class AriaEcContactMap(object):
if not self.settings.contactmap.args.get("onlyreport", False):
self.settings.make_infra()
# ----------------------------- Input -------------------------------- #
if self.settings.contactmap.args.get("prefix"):
self.outprefix = "_".join((get_filename(
self.settings.contactmap.args.get("seq", None)),
self.settings.contactmap.args.get("prefix", "")))
if self.settings.contactmap.args.get("prefix") and len(
self.settings.contactmap.args.get("infiles")) == len(
self.settings.contactmap.args.get("prefix")):
self.outprefix = self.settings.contactmap.args.get("prefix", "")
else:
self.outprefix = get_filename(
self.settings.contactmap.args.get("seq", None))
self.outprefix = get_filename(self.settings.contactmap.args.get("seq",
None))
# Load Sequence file
self.protein.set_aa_sequence(self.settings.contactmap.args.get("seq", None))
# Load secondary structure prediction file
......@@ -70,11 +70,13 @@ class AriaEcContactMap(object):
logger.info("%s map set as reference" % fo.filetype.capitalize())
self.refmap = fo.mapdict
self.reftype = fo.filetype
self.refname = fo.filename
self.refname = fo.filename if type(self.outprefix) != list \
else self.outprefix[idx]
if not self.settings.contactmap.args.get("nofilter"):
self.filter(fo.mapdict, fo.filetype, fo.contactlist,
self.protein, clashlist=fo.clashlist,
outprefix=self.outprefix,
outprefix=self.outprefix[idx] if type(
self.outprefix) == list else self.outprefix,
outdir=self.settings.outdir)
# else:
# Use only position filter
......@@ -82,7 +84,9 @@ class AriaEcContactMap(object):
# self.protein, clashlist=fo.clashlist,
# outprefix=self.outprefix,
# outdir=self.settings.outdir, mapfilters="pos")
self.allresmap[(fo.filename, fo.filetype)] = fo.mapdict
self.allresmap[(fo.filename if type(self.outprefix) != list else
self.outprefix[idx], fo.filetype,
fo.filepath)] = fo.mapdict
try:
refmap = self.refmap["contactmap"]
......@@ -99,7 +103,8 @@ class AriaEcContactMap(object):
if self.settings.contactmap.config.get("save_fig") and not \
self.settings.contactmap.args.get("onlyreport", False):
refmap.saveplot(outdir=outdir,
outprefix="_".join((self.outprefix, "pdb")),
outprefix=self.outprefix[0] if type(
self.outprefix) == list else self.outprefix,
**plotparams)
if self.settings.contactmap.args.get("merge", None):
......@@ -120,20 +125,20 @@ class AriaEcContactMap(object):
self.allresmap[mergekey] = {}
self.allresmap[mergekey]["contactmap"] = up_map
for mapname, mapt in self.allresmap.keys():
for mapname, mapt, mapath in self.allresmap.keys():
prefix = "_".join((self.outprefix, mapname, self.refname))
prefix = "_".join((mapname, self.refname))
if mapname == self.refname:
if not self.settings.contactmap.args.get("onlyreport", False):
refmap.write_contacts(mapname,
outdir=outdir,
scoremap=self.refmap.get("scoremap",
None))
continue
if self.settings.contactmap.args.get("onlyreport", False):
refmap.write_contacts("_".join((self.outprefix, mapt)),
outdir=outdir,
scoremap=self.refmap.get("scoremap",
None))
continue
scoremap = self.allresmap[(mapname, mapt)].get('scoremap', None)
scoremap = self.allresmap[(mapname, mapt, mapath)].get(
'scoremap', None)
# if self.allresmap[mapt].get("contactmap") is not None and \
# self.allresmap[mapt].get("scoremap") is not None:
# Get top contact map/list
......@@ -143,10 +148,9 @@ class AriaEcContactMap(object):
# human_idx=True)[:nb_c]
# elif self.allresmap[mapt].get("contactmap") is not None:
# If no score given, use all contact list
cmpmap = self.allresmap[(mapname, mapt)]["contactmap"]
cmplist = self.allresmap[(mapname, mapt)][
'contactmap'].contact_list(
human_idx=True)
cmpmap = self.allresmap[(mapname, mapt, mapath)]["contactmap"]
cmplist = self.allresmap[(mapname, mapt, mapath)][
'contactmap'].contact_list(human_idx=True)
# else:
# logger.warning("%s map doesn't have any score related. Can't "
# "define top list related to this map" % mapt)
......@@ -165,7 +169,7 @@ class AriaEcContactMap(object):
# TODO: elementwise error with compare method
# Write cmp stats
if not self.settings.contactmap.args.get("onlyreport", False):
cmpmap.write_contacts("_".join((self.outprefix, mapt)),
cmpmap.write_contacts(mapname,
scoremap=scoremap,
outdir=outdir)
cmpmap.compare_contactmap(refmap, cmplist, prefix,
......
No preview for this file type
......@@ -200,11 +200,12 @@ class ProteinMap(Map):
minticks = tickmin(self, shift=1) # Nb graduations
self._maplot = sns.heatmap(self, square=True, cbar=False,
linewidths=1, vmax=1, vmin=-1,
linewidths=0.5, vmax=1, vmin=-1,
cmap=sns.diverging_palette(20, 220, n=7,
as_cmap=True),
xticklabels=minticks[0],
yticklabels=minticks[1])
yticklabels=minticks[1],
linecolor="grey")
return self._maplot
def saveplot(self, outdir='', outprefix="protein", size_fig=10,
......@@ -269,7 +270,7 @@ class ProteinMap(Map):
protmap.__class__.__name__, self.__class__.__name__))
return None
else:
cmplist = protmap.contact_list()
cmplist = protmap.contact_list(human_idx=True)
ymax = len(self.sequence) - 1
if protmap.contact_flags:
......@@ -286,19 +287,23 @@ class ProteinMap(Map):
color = pal[i]
mark = Line2D.filled_markers[i]
for x, y in zip(xind, yind):
self.maplot.axes.scatter(x, y, s=10, c=color,
linewidths=0, alpha=alpha,
self.maplot.axes.scatter(x, y, s=8, c=color,
linewidths=0.1, alpha=alpha,
marker=mark)
else:
logger.info("Contact list: %s" % cmplist)
xind = [x - .5 for x in
zip(*cmplist)[0] + zip(*cmplist)[1]]
yind = [ymax - y + 1.5 for y in
zip(*cmplist)[1] + zip(*cmplist)[0]]
logger.debug("Xind: %s" % xind)
logger.debug("Yind: %s" % yind)
color = "red"
# width = [0.3 for _ in xind]
# for x, y, h in zip(xind, yind, width):
for x, y in zip(xind, yind):
self.maplot.axes.scatter(x, y, s=10, c=color, linewidths=0,
self.maplot.axes.scatter(x, y, s=30, c=color,
linewidths=0,
alpha=alpha)
if save_fig:
self.saveplot(**kwargs)
......@@ -1239,7 +1244,9 @@ class MapFilter:
# Contactmap always filtered
# TODO: could set a treshold instead of n_factor
nb_c = int(len(mapdict["contactmap"].sequence) * int(
logger.info("Setting contact number with treshold %s" %
self.settings.get("n_factor"))
nb_c = int(len(mapdict["contactmap"].sequence) * float(
self.settings.get("n_factor")))
nb_c = nb_c if nb_c < len(mapdict["contactmap"].contactset()) else len(
mapdict["contactmap"].contactset())
......@@ -1292,8 +1299,8 @@ class MapFilter:
for clash_t in sorted(clash_dict.keys()):
if clash_dict[clash_t] and raw_contact in clash_dict[clash_t]:
clash = meta_clash[clash_t]["flag"]
meta_clash[clash_t]["msg"] += """\
{clash_type} flag at pair {pair_nb} : res {res1} and res {res2} {clash_desc}
meta_clash[clash_t]["msg"] += """
{clash_type} flag at pair {pair_nb} : res {res1} and res {res2} {clash_desc}\
""".format(clash_type=clash, pair_nb=icontact + 1,
clash_desc=desc_dict.get(raw_contact, ''),
res1=contact[0], res2=contact[1])
......@@ -1305,9 +1312,9 @@ class MapFilter:
else:
op = "added"
ctype = clash
meta_clash[clash_t]["warn"] += r"\n/!\ Clash: {clash_desc} {" \
r"clash} flag for contact " \
r"{res_pos} ({res1}, {res2})".format(
meta_clash[clash_t]["warn"] += "\n/!\ Clash: {clash_desc} {" \
"clash} flag for contact " \
"{res_pos} ({res1}, {res2})".format(
clash_desc=op, clash=ctype, res_pos=icontact + 1,
res1=contact[0], res2=contact[1])
......@@ -1315,6 +1322,7 @@ class MapFilter:
titleprint(out, progname=__doc__, desc='Contacts filter')
for flt in ("nd", "pos", "cons", "ssclash", "cys"):
out.write('''
{filter_desc}
{hd}
'''.format(filter_desc=meta_clash[flt]["desc"].capitalize(),
......
No preview for this file type
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment