Coverage for biobb_dna/stiffness/basepair_stiffness.py: 78%

1#!/usr/bin/env python3 

2 

3"""Module containing the HelParStiffness class and the command line interface.""" 

4 

5import argparse 

6from pathlib import Path 

7from typing import Optional 

8 

9import matplotlib.pyplot as plt 

10import numpy as np 

11import pandas as pd 

12from biobb_common.configuration import settings 

13from biobb_common.generic.biobb_object import BiobbObject 

14from biobb_common.tools.file_utils import launchlogger 

15 

16from biobb_dna.utils.common import _from_string_to_list 

17from biobb_dna.utils.loader import load_data 

18 

19 

20class BPStiffness(BiobbObject): 

21 """ 

22 | biobb_dna BPStiffness 

23 | Calculate stiffness constants matrix between all six helical parameters for a single base pair step. 

24 | Calculate stiffness constants matrix between all six helical parameters for a single base pair step. 

25 

26 Args: 

27 input_filename_shift (str): Path to csv file with data for helical parameter 'shift'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/stiffness/series_shift_AA.csv>`_. Accepted formats: csv (edam:format_3752) 

28 input_filename_slide (str): Path to csv file with data for helical parameter 'slide'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/stiffness/series_slide_AA.csv>`_. Accepted formats: csv (edam:format_3752) 

29 input_filename_rise (str): Path to csv file with data for helical parameter 'rise'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/stiffness/series_rise_AA.csv>`_. Accepted formats: csv (edam:format_3752) 

30 input_filename_tilt (str): Path to csv file with data for helical parameter 'tilt'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/stiffness/series_tilt_AA.csv>`_. Accepted formats: csv (edam:format_3752) 

31 input_filename_roll (str): Path to csv file with data for helical parameter 'roll'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/stiffness/series_roll_AA.csv>`_. Accepted formats: csv (edam:format_3752) 

32 input_filename_twist (str): Path to csv file with data for helical parameter 'twist'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/stiffness/series_twist_AA.csv>`_. Accepted formats: csv (edam:format_3752) 

33 output_csv_path (str): Path to directory where stiffness matrix file is saved as a csv file. File type: output. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/reference/stiffness/stiffbp_ref.csv>`_. Accepted formats: csv (edam:format_3752) 

34 output_jpg_path (str): Path to directory where stiffness heatmap image is saved as a jpg file. File type: output. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/reference/stiffness/stiffbp_ref.jpg>`_. Accepted formats: jpg (edam:format_3579) 

35 properties (dict): 

36 * **KT** (*float*) - (0.592186827) Value of Boltzmann temperature factor. 

37 * **scaling** (*list*) - ([1, 1, 1, 10.6, 10.6, 10.6]) Values by which to scale stiffness. Positions correspond to helical parameters in the order: shift, slide, rise, tilt, roll, twist. 

38 * **remove_tmp** (*bool*) - (True) [WF property] Remove temporal files. 

39 * **restart** (*bool*) - (False) [WF property] Do not execute if output files exist. 

40 * **sandbox_path** (*str*) - ("./") [WF property] Parent path to the sandbox directory. 

41 

42 Examples: 

43 This is a use example of how to use the building block from Python:: 

44 

45 from biobb_dna.stiffness.basepair_stiffness import basepair_stiffness 

46 

47 prop = { 

48 'KT': 0.592186827, 

49 'scaling': [1, 1, 1, 10.6, 10.6, 10.6] 

50 } 

51 basepair_stiffness( 

52 input_filename_shift='path/to/basepair/shift.csv', 

53 input_filename_slide='path/to/basepair/slide.csv', 

54 input_filename_rise='path/to/basepair/rise.csv', 

55 input_filename_tilt='path/to/basepair/tilt.csv', 

56 input_filename_roll='path/to/basepair/roll.csv', 

57 input_filename_twist='path/to/basepair/twist.csv', 

58 output_csv_path='path/to/output/file.csv', 

59 output_jpg_path='path/to/output/plot.jpg', 

60 properties=prop) 

61 Info: 

62 * wrapped_software: 

63 * name: In house 

64 * license: Apache-2.0 

65 * ontology: 

66 * name: EDAM 

67 * schema: http://edamontology.org/EDAM.owl 

68 

69 """ 

70 

71 def __init__( 

72 self, 

73 input_filename_shift, 

74 input_filename_slide, 

75 input_filename_rise, 

76 input_filename_tilt, 

77 input_filename_roll, 

78 input_filename_twist, 

79 output_csv_path, 

80 output_jpg_path, 

81 properties=None, 

82 **kwargs, 

83 ) -> None: 

84 properties = properties or {} 

85 

86 # Call parent class constructor 

87 super().__init__(properties) 

88 self.locals_var_dict = locals().copy() 

89 

90 # Input/Output files 

91 self.io_dict = { 

92 "in": { 

93 "input_filename_shift": input_filename_shift, 

94 "input_filename_slide": input_filename_slide, 

95 "input_filename_rise": input_filename_rise, 

96 "input_filename_tilt": input_filename_tilt, 

97 "input_filename_roll": input_filename_roll, 

98 "input_filename_twist": input_filename_twist, 

99 }, 

100 "out": { 

101 "output_csv_path": output_csv_path, 

102 "output_jpg_path": output_jpg_path, 

103 }, 

104 } 

105 

106 self.properties = properties 

107 self.KT = properties.get("KT", 0.592186827) 

108 self.scaling = [ 

109 int(elem) 

110 for elem in _from_string_to_list( 

111 properties.get("scaling", [1, 1, 1, 10.6, 10.6, 10.6]) 

112 ) 

113 ] 

114 

115 # Check the properties 

116 self.check_properties(properties) 

117 self.check_arguments() 

118 

119 @launchlogger 

120 def launch(self) -> int: 

121 """Execute the :class:`BPStiffness <stiffness.basepair_stiffness.BPStiffness>` object.""" 

122 

123 # Setup Biobb 

124 if self.check_restart(): 

125 return 0 

126 self.stage_files() 

127 

128 # read input 

129 shift = load_data(self.stage_io_dict["in"]["input_filename_shift"]) 

130 slide = load_data(self.stage_io_dict["in"]["input_filename_slide"]) 

131 rise = load_data(self.stage_io_dict["in"]["input_filename_rise"]) 

132 tilt = load_data(self.stage_io_dict["in"]["input_filename_tilt"]) 

133 roll = load_data(self.stage_io_dict["in"]["input_filename_roll"]) 

134 twist = load_data(self.stage_io_dict["in"]["input_filename_twist"]) 

135 

136 # build matrix cols_arr from helpar input data files 

137 coordinates = ["shift", "slide", "rise", "tilt", "roll", "twist"] 

138 basepairname = shift.columns[0] 

139 helpar_matrix = pd.concat([shift, slide, rise, tilt, roll, twist], axis=1) 

140 helpar_matrix.columns = coordinates 

141 # covariance 

142 cov_df = helpar_matrix.cov() 

143 # stiffness 

144 stiff = np.linalg.inv(cov_df) * self.KT 

145 stiff_diag = stiff * np.array(self.scaling) 

146 stiff_df = pd.DataFrame(stiff_diag, columns=cov_df.columns, index=cov_df.index) 

147 stiff_df.index.name = basepairname 

148 

149 # save csv data 

150 stiff_df.to_csv(Path(self.stage_io_dict["out"]["output_csv_path"])) 

151 

152 # create heatmap 

153 fig, axs = plt.subplots(1, 1, dpi=300, tight_layout=True) 

154 axs.pcolor(stiff_df) 

155 # Loop over data dimensions and create text annotations. 

156 for i in range(len(stiff_df)): 

157 for j in range(len(stiff_df)): 

158 axs.text( 

159 j + 0.5, 

160 i + 0.5, 

161 f"{stiff_df[coordinates[j]].loc[coordinates[i]]:.2f}", 

162 ha="center", 

163 va="center", 

164 color="w", 

165 ) 

166 axs.text( 

167 0, 

168 -1.35, 

169 "Units:\n" 

170 "Diagonal Shift/Slide/Rise in kcal/(mol*Ų), Diagonal Tilt/Roll/Twist in kcal/(mol*degree²)\n" 

171 "Out of Diagonal: Shift/Slide/Rise in kcal/(mol*Å), Out of Diagonal Tilt/Roll/Twist in kcal/(mol*degree)", 

172 fontsize=6, 

173 ) 

174 axs.set_xticks([i + 0.5 for i in range(len(stiff_df))]) 

175 axs.set_xticklabels(stiff_df.columns, rotation=90) 

176 axs.set_yticks([i + 0.5 for i in range(len(stiff_df))]) 

177 axs.set_yticklabels(stiff_df.index) 

178 axs.set_title(f"Stiffness Constants for Base Pair Step '{basepairname}'") 

179 fig.tight_layout() 

180 fig.savefig(self.stage_io_dict["out"]["output_jpg_path"], format="jpg") 

181 plt.close() 

182 

183 # Copy files to host 

184 self.copy_to_host() 

185 

186 # Remove temporary file(s) 

187 # self.tmp_files.extend([self.stage_io_dict.get("unique_dir", "")]) 

188 self.remove_tmp_files() 

189 

190 self.check_arguments(output_files_created=True, raise_exception=False) 

191 

192 return 0 

193 

194 

195def basepair_stiffness( 

196 input_filename_shift: str, 

197 input_filename_slide: str, 

198 input_filename_rise: str, 

199 input_filename_tilt: str, 

200 input_filename_roll: str, 

201 input_filename_twist: str, 

202 output_csv_path: str, 

203 output_jpg_path: str, 

204 properties: Optional[dict] = None, 

205 **kwargs, 

206) -> int: 

207 """Create :class:`BPStiffness <stiffness.basepair_stiffness.BPStiffness>` class and 

208 execute the :meth:`launch() <stiffness.basepair_stiffness.BPStiffness.BPStiffness.launch>` method.""" 

209 

210 return BPStiffness( 

211 input_filename_shift=input_filename_shift, 

212 input_filename_slide=input_filename_slide, 

213 input_filename_rise=input_filename_rise, 

214 input_filename_tilt=input_filename_tilt, 

215 input_filename_roll=input_filename_roll, 

216 input_filename_twist=input_filename_twist, 

217 output_csv_path=output_csv_path, 

218 output_jpg_path=output_jpg_path, 

219 properties=properties, 

220 **kwargs, 

221 ).launch() 

222 

223 basepair_stiffness.__doc__ = BPStiffness.__doc__ 

224 

225 

226def main(): 

227 """Command line execution of this building block. Please check the command line documentation.""" 

228 parser = argparse.ArgumentParser( 

229 description="Calculate stiffness constants matrix between all six helical parameters for a single base pair step.", 

230 formatter_class=lambda prog: argparse.RawTextHelpFormatter(prog, width=99999), 

231 ) 

232 parser.add_argument("--config", required=False, help="Configuration file") 

233 

234 required_args = parser.add_argument_group("required arguments") 

235 required_args.add_argument( 

236 "--input_filename_shift", 

237 required=True, 

238 help="Path to csv file with shift inputs. Accepted formats: csv.", 

239 ) 

240 required_args.add_argument( 

241 "--input_filename_slide", 

242 required=True, 

243 help="Path to csv file with slide inputs. Accepted formats: csv.", 

244 ) 

245 required_args.add_argument( 

246 "--input_filename_rise", 

247 required=True, 

248 help="Path to csv file with rise inputs. Accepted formats: csv.", 

249 ) 

250 required_args.add_argument( 

251 "--input_filename_tilt", 

252 required=True, 

253 help="Path to csv file with tilt inputs. Accepted formats: csv.", 

254 ) 

255 required_args.add_argument( 

256 "--input_filename_roll", 

257 required=True, 

258 help="Path to csv file with roll inputs. Accepted formats: csv.", 

259 ) 

260 required_args.add_argument( 

261 "--input_filename_twist", 

262 required=True, 

263 help="Path to csv file with twist inputs. Accepted formats: csv.", 

264 ) 

265 required_args.add_argument( 

266 "--output_csv_path", 

267 required=True, 

268 help="Path to output covariance data file. Accepted formats: csv.", 

269 ) 

270 required_args.add_argument( 

271 "--output_jpg_path", 

272 required=True, 

273 help="Path to output covariance data file. Accepted formats: csv.", 

274 ) 

275 

276 args = parser.parse_args() 

277 args.config = args.config or "{}" 

278 properties = settings.ConfReader(config=args.config).get_prop_dic() 

279 

280 basepair_stiffness( 

281 input_filename_shift=args.input_filename_shift, 

282 input_filename_slide=args.input_filename_slide, 

283 input_filename_rise=args.input_filename_rise, 

284 input_filename_tilt=args.input_filename_tilt, 

285 input_filename_roll=args.input_filename_roll, 

286 input_filename_twist=args.input_filename_twist, 

287 output_csv_path=args.output_csv_path, 

288 output_jpg_path=args.output_jpg_path, 

289 properties=properties, 

290 ) 

291 

292 

293if __name__ == "__main__": 

294 main()