Mercurial > repos > chemteam > mdanalysis_extract_rmsd

comparison pca_cosine.py @ 0:743bd6aa3c7a draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload for repository https://github.com/galaxycomputationalchemistry/galaxy-tools-compchem/ commit 1b23e024af45cc0999d9142d07de6897d4189ec2"
author chemteam
date Mon, 24 Aug 2020 16:41:41 -0400
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:743bd6aa3c7a
1 #!/usr/bin/env python
2
3 import argparse
4 import csv
5 import sys
6
7 import MDAnalysis as mda
8 import MDAnalysis.analysis.pca as pca
9
10 import numpy as np
11
12
13 def parse_command_line(argv):
14 parser = argparse.ArgumentParser()
15 parser.add_argument('--itraj', help='input traj')
16 parser.add_argument('--istr', help='input str')
17 parser.add_argument('--itrajext', help='input traj ext')
18 parser.add_argument('--istrext', help='input str ext')
19 parser.add_argument('--icomponents', help='number of principle components')
20 parser.add_argument('--iindex', help='index of the PC')
21 parser.add_argument('--output', help='output')
22 parser.add_argument('--cosout', help='cosine output')
23 return parser.parse_args()
24
25
26 args = parse_command_line(sys.argv)
27
28 u = mda.Universe(args.istr, args.itraj,
29 topology_format=args.istrext, format=args.itrajext)
30
31 components = int(args.icomponents)
32 pca_index = int(args.iindex)
33
34 PSF_pca = pca.PCA(u, select='backbone')
35 PSF_pca.run()
36 n_pcs = np.where(PSF_pca.cumulated_variance > 0.95)[0][0]
37 atomgroup = u.select_atoms('backbone')
38
39 pca_space = PSF_pca.transform(atomgroup, n_components=components)
40 cosine = mda.analysis.pca.cosine_content(pca_space, pca_index)
41
42 PCA = list(pca_space)
43
44 with open(args.output, 'w') as f:
45 writer = csv.writer(f, delimiter='\t')
46 writer.writerows(PCA)
47
48 with open(args.cosout, 'w') as f1:
49 f1.write(str(cosine))