Coverage for mddb_workflow / analyses / rgyr.py: 93%

27 statements  

« prev     ^ index     » next       coverage.py v7.12.0, created at 2025-12-03 18:45 +0000

1from os import remove 

2from numpy import mean, std 

3 

4from mddb_workflow.tools.get_reduced_trajectory import get_reduced_trajectory 

5from mddb_workflow.tools.xvg_parse import xvg_parse 

6from mddb_workflow.utils.auxiliar import save_json 

7from mddb_workflow.utils.constants import OUTPUT_RGYR_FILENAME 

8from mddb_workflow.utils.gmx_spells import run_gromacs 

9from mddb_workflow.utils.type_hints import * 

10 

11# Set an auxiliar data filename 

12rgyr_data_filename = '.rgyr_data.xvg' 

13 

14 

15def rgyr( 

16 structure_file: 'File', 

17 trajectory_file: 'File', 

18 output_directory: str, 

19 snapshots: int, 

20 structure: 'Structure', 

21 pbc_selection: 'Selection', 

22 frames_limit: int = 5000 

23): 

24 """Perform the RMSd analysis. Use the first trajectory frame in .pdb format as a reference.""" 

25 # Set the main output filepath 

26 output_analysis_filepath = f'{output_directory}/{OUTPUT_RGYR_FILENAME}' 

27 

28 # Use a reduced trajectory in case the original trajectory has many frames 

29 reduced_trajectory_filepath, step, frames = get_reduced_trajectory( 

30 structure_file, 

31 trajectory_file, 

32 snapshots, 

33 frames_limit, 

34 ) 

35 

36 # Generate a custom index file to exclude PBC residues from the analysis 

37 not_pbc_selection = structure.invert_selection(pbc_selection) 

38 if not not_pbc_selection: 

39 print(' No selection to run the analysis') 

40 return 

41 selection_name = 'pbc_atoms' 

42 ndx_selection = not_pbc_selection.to_ndx(selection_name) 

43 ndx_filename = '.not_pbc.ndx' 

44 with open(ndx_filename, 'w') as file: 

45 file.write(ndx_selection) 

46 

47 # Run Gromacs 

48 run_gromacs(f'gyrate -s {structure_file.path} -f {reduced_trajectory_filepath} \ 

49 -o {rgyr_data_filename} -n {ndx_filename} -mode geometry', user_input=selection_name) 

50 

51 # Read the output file and parse it 

52 raw_rgyr_data = xvg_parse(rgyr_data_filename, ['times', 'rgyr', 'rgyrx', 'rgyry', 'rgyrz']) 

53 

54 # Format data 

55 rgyr_data = { 

56 'start': 0, 

57 'step': step, 

58 'y': { 

59 'rgyr': { 

60 'average': mean(raw_rgyr_data['rgyr']), 

61 'stddev': std(raw_rgyr_data['rgyr']), 

62 'min': min(raw_rgyr_data['rgyr']), 

63 'max': max(raw_rgyr_data['rgyr']), 

64 'data': raw_rgyr_data['rgyr'] 

65 }, 

66 'rgyrx': { 

67 'average': mean(raw_rgyr_data['rgyrx']), 

68 'stddev': std(raw_rgyr_data['rgyrx']), 

69 'min': min(raw_rgyr_data['rgyrx']), 

70 'max': max(raw_rgyr_data['rgyrx']), 

71 'data': raw_rgyr_data['rgyrx'] 

72 }, 

73 'rgyry': { 

74 'average': mean(raw_rgyr_data['rgyry']), 

75 'stddev': std(raw_rgyr_data['rgyry']), 

76 'min': min(raw_rgyr_data['rgyry']), 

77 'max': max(raw_rgyr_data['rgyry']), 

78 'data': raw_rgyr_data['rgyry'] 

79 }, 

80 'rgyrz': { 

81 'average': mean(raw_rgyr_data['rgyrz']), 

82 'stddev': std(raw_rgyr_data['rgyrz']), 

83 'min': min(raw_rgyr_data['rgyrz']), 

84 'max': max(raw_rgyr_data['rgyrz']), 

85 'data': raw_rgyr_data['rgyrz'] 

86 } 

87 }, 

88 'version': '0.0.1', 

89 } 

90 

91 # Export formatted data to a json file 

92 save_json(rgyr_data, output_analysis_filepath) 

93 

94 # Remove residual files 

95 remove(rgyr_data_filename) 

96 remove(ndx_filename)