3 from __future__
import print_function
5 from IMP
import ArgumentParser
7 __doc__ =
"Score each of a set of combinations."
14 colors[
"Rpt1"] = [0.78, 0.78, 0.73]
15 colors[
"Rpt2"] = [0.78, 0.66, 0.58]
16 colors[
"Rpt3"] = [0.77, 0.43, 0.5]
17 colors[
"Rpt4"] = [0.76, 0.29, 0.67]
18 colors[
"Rpt5"] = [0.51, 0.14, 0.75]
19 colors[
"Rpt6"] = [0.0, 0., 0.75]
20 colors[
"Rpn1"] = [0.34, 0.36, 0.27]
21 colors[
"Rpn2"] = [0.42, 0.43, 0.36]
22 colors[
"Rpn3"] = [0.49, 0.5, 0.44]
23 colors[
"Rpn5"] = [0.56, 0.57, 0.51]
24 colors[
"Rpn6"] = [0.64, 0.64, 0.59]
25 colors[
"Rpn7"] = [0.71, 0.71, 0.66]
26 colors[
"Rpn8"] = [0.78, 0.78, 0.74]
27 colors[
"Rpn9"] = [1, 0, 0]
28 colors[
"Rpn10"] = [0, 1, 0]
29 colors[
"Rpn11"] = [0, 0, 1]
30 colors[
"Rpn12"] = [0.5, 0.2, 0.4]
31 colors[
"a1"] = [0.78, 0.78, 0.73]
32 colors[
"a2"] = [0.78, 0.66, 0.58]
33 colors[
"a3"] = [0.77, 0.43, 0.5]
34 colors[
"a4"] = [0.76, 0.29, 0.67]
35 colors[
"a5"] = [0.51, 0.14, 0.75]
36 colors[
"a6"] = [0.0, 0., 0.75]
37 colors[
"a7"] = [0.34, 0.36, 0.27]
38 colors[
"a8"] = [0.42, 0.43, 0.36]
39 colors[
"a9"] = [0.49, 0.5, 0.44]
40 colors[
"a10"] = [0.56, 0.57, 0.51]
42 colors[
"a11"] = [0.78, 0.78, 0.73]
43 colors[
"a12"] = [0.78, 0.66, 0.58]
44 colors[
"a13"] = [0.77, 0.43, 0.5]
45 colors[
"a14"] = [0.76, 0.29, 0.67]
46 colors[
"a15"] = [0.51, 0.14, 0.75]
47 colors[
"a16"] = [0.0, 0., 0.75]
48 colors[
"a17"] = [0.34, 0.36, 0.27]
49 colors[
"a18"] = [0.42, 0.43, 0.36]
50 colors[
"a19"] = [0.49, 0.5, 0.44]
51 colors[
"a20"] = [0.56, 0.57, 0.51]
53 colors[
"a21"] = [0.78, 0.78, 0.73]
54 colors[
"a22"] = [0.78, 0.66, 0.58]
55 colors[
"a23"] = [0.77, 0.43, 0.5]
56 colors[
"a24"] = [0.76, 0.29, 0.67]
57 colors[
"a25"] = [0.51, 0.14, 0.75]
58 colors[
"a26"] = [0.0, 0., 0.75]
59 colors[
"a27"] = [0.34, 0.36, 0.27]
60 colors[
"a28"] = [0.42, 0.43, 0.36]
61 colors[
"a29"] = [0.49, 0.5, 0.44]
62 colors[
"a30"] = [0.56, 0.57, 0.51]
66 def decompose(dmap, mhs):
71 full_sampled_map.set_particles(all_ps)
72 full_sampled_map.resample()
73 full_sampled_map.calcRMS()
75 dmap.get_number_of_voxels(
77 ).dmean * full_sampled_map.get_header(
80 lower = dmap.get_number_of_voxels(
82 ).rms * full_sampled_map.get_header(
84 norm_factors = [upper, lower]
85 print(
"===============my norm factors:", upper, lower)
89 def score_each_protein(dmap, mhs, sd):
90 norm_factors = decompose(dmap, mhs)
93 for i
in range(len(mhs)):
97 mh_dmap.set_particles(leaves)
101 sd.get_component_header(i).get_transformations_fn())
103 for fit
in fits[:15]:
114 scores.append(mh_scores)
115 print(
"=====mol", i, mh_scores)
120 desc =
"Score each of a set of combinations."
121 p = ArgumentParser(description=desc)
122 p.add_argument(
"-m",
"--max", dest=
"max", type=int, default=999999999,
123 help=
"maximum number of fits considered")
124 p.add_argument(
"assembly_file", help=
"assembly file name")
125 p.add_argument(
"proteomics_file", help=
"proteomics file name")
126 p.add_argument(
"mapping_file", help=
"mapping file name")
127 p.add_argument(
"param_file", help=
"parameter file name")
128 p.add_argument(
"combinations_file", help=
"combinations file name")
129 p.add_argument(
"scores_file", help=
"output scores file name")
130 return p.parse_args()
133 def run(asmb_fn, proteomics_fn, mapping_fn, params_fn, combs_fn,
134 scored_comb_output_fn, max_comb):
136 dmap = IMP.em.read_map(asmb.get_assembly_header().get_dens_fn())
137 dmap.get_header().set_resolution(
139 dmap.update_voxel_size(asmb.get_assembly_header().get_spacing())
140 dmap.set_origin(asmb.get_assembly_header().get_origin())
141 threshold = asmb.get_assembly_header().get_threshold()
147 alignment_params = IMP.multifit.AlignmentParams(params_fn)
148 alignment_params.show()
151 print(
"=========", combs_fn)
159 prot_data, mapping_fn)
161 _ = mapping_data.get_anchors()
167 mapping_data, asmb, alignment_params)
168 align.set_fast_scoring(
False)
170 mdl = align.get_model()
171 mhs = align.get_molecules()
172 align.add_states_and_filters()
176 align.set_density_map(dmap, threshold)
177 for i, mh
in enumerate(mhs):
178 ensmb.add_component_and_fits(
180 asmb.get_component_header(i).get_transformations_fn()))
183 rgb = colors[mh.get_name()]
185 rgb = colors[names[i]]
186 color = IMP.display.Color(rgb[0], rgb[1], rgb[2])
187 for p in IMP.core.get_leaves(mh):
188 g= IMP.display.XYZRGeometry(p)
196 align.add_all_restraints()
198 rs = align.get_restraint_set().get_restraints()
199 print(
"Get number of restraints:", len(rs))
201 rr = IMP.RestraintSet.get_from(r)
202 for i
in range(rr.get_number_of_restraints()):
203 print(rr.get_restraint(i).get_name())
204 output = open(scored_comb_output_fn,
"w")
208 for i
in range(asmb.get_number_of_component_headers()):
209 c = asmb.get_component_header(i)
210 fn = c.get_reference_fn()
215 rr = IMP.RestraintSet.get_from(r)
216 for i
in range(rr.get_number_of_restraints()):
217 output.write(rr.get_restraint(i).get_name() +
"|")
222 print(
"Number of combinations:", len(combs[:max_comb]))
224 print(
"native score")
227 rr = IMP.RestraintSet.get_from(r)
228 for j
in range(rr.get_number_of_restraints()):
229 print(rr.get_restraint(j).get_name(), rr.evaluate(
False))
232 for i, comb
in enumerate(combs[:max_comb]):
233 print(
"Scoring combination:", comb)
234 ensmb.load_combination(comb)
237 rr = IMP.RestraintSet.get_from(r)
238 for j
in range(rr.get_number_of_restraints()):
239 print(rr.get_restraint(j).get_name())
240 rscore = rr.evaluate(
False)
242 num_violated = num_violated + 1
244 print(str(all_leaves[0]) +
" :: " + str(all_leaves[-1]))
245 score = sf.evaluate(
False)
247 msg =
"COMB" + str(i) +
"|"
249 rr = IMP.RestraintSet.get_from(r)
250 for j
in range(rr.get_number_of_restraints()):
251 current_name = rr.get_restraint(j).get_name()
252 if current_name != prev_name:
253 msg +=
' ' + current_name +
' '
254 prev_name = current_name
255 rscore = rr.get_restraint(j).
evaluate(
False)
256 msg += str(rscore) +
"|"
258 num_violated = num_violated + 1
262 num_violated) +
"||||" + str(
263 fitr.evaluate(
False)) +
"||:"
267 output.write(msg +
"\n")
269 ensmb.unload_combination(comb)
275 run(args.assembly_file, args.proteomics_file, args.mapping_file,
276 args.param_file, args.combinations_file, args.scores_file, args.max)
279 if __name__ ==
"__main__":
An ensemble of fitting solutions.
double get_coarse_cc_coefficient(const DensityMap *grid1, const DensityMap *grid2, double grid2_voxel_data_threshold, bool allow_padding=false, FloatPair norm_factors=FloatPair(0., 0.))
Calculates the cross correlation coefficient between two maps.
void write_pdb(const Selection &mhd, TextOutput out, unsigned int model=1)
def evaluate
Evaluate the score of the restraint.
Create a scoring function on a list of restraints.
SettingsData * read_settings(const char *filename)
GenericHierarchies get_leaves(Hierarchy mhd)
Get all the leaves of the bit of hierarchy.
void read_pdb(TextInput input, int model, Hierarchy h)
ProteinsAnchorsSamplingSpace read_protein_anchors_mapping(multifit::ProteomicsData *prots, const std::string &anchors_prot_map_fn, int max_paths=INT_MAX)
Align proteomics graph to EM density map.
Class for sampling a density map from particles.
double get_rmsd(const Selection &s0, const Selection &s1)
void transform(XYZ a, const algebra::Transformation3D &tr)
Apply a transformation to the particle.
Fitting atomic structures into a cryo-electron microscopy density map.
ProteomicsData * read_proteomics_data(const char *proteomics_fn)
Proteomics reader.
void set_log_level(LogLevel l)
Set the current global log level.
IntsList read_paths(const char *txt_filename, int max_paths=INT_MAX)
Read paths.
Calculate score based on fit to EM map.
FittingSolutionRecords read_fitting_solutions(const char *fitting_fn)
Fitting solutions reader.
double get_resolution(Model *m, ParticleIndex pi)
Estimate the resolution of the hierarchy as used by Representation.