Coverage for biobb_dna / intrabp_correlations / intrahpcorr.py: 91%

1#!/usr/bin/env python3 

2 

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

4from typing import Optional 

5 

6import pandas as pd 

7import numpy as np 

8import matplotlib.pyplot as plt 

9 

10from biobb_common.generic.biobb_object import BiobbObject 

11from biobb_common.tools.file_utils import launchlogger 

12from biobb_dna.utils.loader import load_data 

13 

14 

15class IntraHelParCorrelation(BiobbObject): 

16 """ 

17 | biobb_dna IntraHelParCorrelation 

18 | Calculate correlation between helical parameters for a single intra-base pair. 

19 | Calculate correlation between helical parameters for a single intra-base pair. 

20 

21 Args: 

22 input_filename_shear (str): Path to .csv file with data for helical parameter 'shear'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/correlation/series_shear_A.csv>`_. Accepted formats: csv (edam:format_3752). 

23 input_filename_stretch (str): Path to .csv file with data for helical parameter 'stretch'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/correlation/series_stretch_A.csv>`_. Accepted formats: csv (edam:format_3752). 

24 input_filename_stagger (str): Path to .csv file with data for helical parameter 'stagger'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/correlation/series_stagger_A.csv>`_. Accepted formats: csv (edam:format_3752). 

25 input_filename_buckle (str): Path to .csv file with data for helical parameter 'buckle'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/correlation/series_buckle_A.csv>`_. Accepted formats: csv (edam:format_3752). 

26 input_filename_propel (str): Path to .csv file with data for helical parameter 'propeller'. File type: input. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/data/correlation/series_propel_A.csv>`_. Accepted formats: csv (edam:format_3752). 

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

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

29 output_jpg_path (str): Path to .jpg file where output is saved. File type: output. `Sample file <https://raw.githubusercontent.com/bioexcel/biobb_dna/master/biobb_dna/test/reference/correlation/intra_hpcorr_ref.jpg>`_. Accepted formats: jpg (edam:format_3579). 

30 properties (dict): 

31 * **base** (*str*) - (None) Name of base analyzed. 

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

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

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

35 

36 Examples: 

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

38 

39 from biobb_dna.intrabp_correlations.intrahpcorr import intrahpcorr 

40 

41 prop = { 

42 'base': 'A', 

43 } 

44 intrahpcorr( 

45 input_filename_shear='path/to/shear.csv', 

46 input_filename_stretch='path/to/stretch.csv', 

47 input_filename_stagger='path/to/stagger.csv', 

48 input_filename_buckle='path/to/buckle.csv', 

49 input_filename_propel='path/to/propel.csv', 

50 input_filename_opening='path/to/opening.csv', 

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

52 output_jpg_path='path/to/output/file.jpg', 

53 properties=prop) 

54 Info: 

55 * wrapped_software: 

56 * name: In house 

57 * license: Apache-2.0 

58 * ontology: 

59 * name: EDAM 

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

61 

62 """ 

63 

64 def __init__( 

65 self, input_filename_shear, input_filename_stretch, 

66 input_filename_stagger, input_filename_buckle, 

67 input_filename_propel, input_filename_opening, 

68 output_csv_path, output_jpg_path, 

69 properties=None, **kwargs) -> None: 

70 properties = properties or {} 

71 

72 # Call parent class constructor 

73 super().__init__(properties) 

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

75 

76 # Input/Output files 

77 self.io_dict = { 

78 'in': { 

79 'input_filename_shear': input_filename_shear, 

80 'input_filename_stretch': input_filename_stretch, 

81 'input_filename_stagger': input_filename_stagger, 

82 'input_filename_buckle': input_filename_buckle, 

83 'input_filename_propel': input_filename_propel, 

84 'input_filename_opening': input_filename_opening 

85 }, 

86 'out': { 

87 'output_csv_path': output_csv_path, 

88 'output_jpg_path': output_jpg_path 

89 } 

90 } 

91 

92 self.properties = properties 

93 self.base = properties.get("base", None) 

94 

95 # Check the properties 

96 self.check_properties(properties) 

97 self.check_arguments() 

98 

99 @launchlogger 

100 def launch(self) -> int: 

101 """Execute the :class:`IntraHelParCorrelation <intrabp_correlations.intrahpcorr.IntraHelParCorrelation>` object.""" 

102 

103 # Setup Biobb 

104 if self.check_restart(): 

105 return 0 

106 self.stage_files() 

107 

108 # read input 

109 shear = load_data(self.stage_io_dict["in"]["input_filename_shear"]) 

110 stretch = load_data(self.stage_io_dict["in"]["input_filename_stretch"]) 

111 stagger = load_data(self.stage_io_dict["in"]["input_filename_stagger"]) 

112 buckle = load_data(self.stage_io_dict["in"]["input_filename_buckle"]) 

113 propel = load_data(self.stage_io_dict["in"]["input_filename_propel"]) 

114 opening = load_data(self.stage_io_dict["in"]["input_filename_opening"]) 

115 

116 # get base 

117 if self.base is None: 

118 self.base = shear.columns[0] 

119 

120 # make matrix 

121 # coordinates = ["shear", "stretch", "stagger", "buckle", "propel", "opening"] 

122 coordinates = [ 

123 "shear", "stretch", "stagger", "buckle", "propel", "opening"] 

124 corr_matrix = pd.DataFrame( 

125 np.eye(6, 6), index=coordinates, columns=coordinates) 

126 

127 # shear 

128 # corr_matrix["shear"]["stretch"] = shear.corrwith(stretch, method="pearson") 

129 corr_matrix.loc["stretch", "shear"] = shear.corrwith(stretch, method="pearson").values[0] 

130 # corr_matrix["shear"]["stagger"] = shear.corrwith(stagger, method="pearson") 

131 corr_matrix.loc["stagger", "shear"] = shear.corrwith(stagger, method="pearson").values[0] 

132 # corr_matrix["shear"]["buckle"] = shear.corrwith(buckle, method=self.circlineal) 

133 corr_matrix.loc["buckle", "shear"] = shear.corrwith(buckle, method=self.circlineal).values[0] # type: ignore 

134 # corr_matrix["shear"]["propel"] = shear.corrwith(propel, method=self.circlineal) 

135 corr_matrix.loc["propel", "shear"] = shear.corrwith(propel, method=self.circlineal).values[0] # type: ignore 

136 # corr_matrix["shear"]["opening"] = shear.corrwith(opening, method=self.circlineal) 

137 corr_matrix.loc["opening", "shear"] = shear.corrwith(opening, method=self.circlineal).values[0] # type: ignore 

138 # symmetric values 

139 # corr_matrix["stretch"]["shear"] = corr_matrix["shear"]["stretch"] 

140 corr_matrix.loc["shear", "stretch"] = corr_matrix.loc["stretch", "shear"] 

141 # corr_matrix["stagger"]["shear"] = corr_matrix["shear"]["stagger"] 

142 corr_matrix.loc["shear", "stagger"] = corr_matrix.loc["stagger", "shear"] 

143 # corr_matrix["buckle"]["shear"] = corr_matrix["shear"]["buckle"] 

144 corr_matrix.loc["shear", "buckle"] = corr_matrix.loc["buckle", "shear"] 

145 # corr_matrix["propel"]["shear"] = corr_matrix["shear"]["propel"] 

146 corr_matrix.loc["shear", "propel"] = corr_matrix.loc["propel", "shear"] 

147 # corr_matrix["opening"]["shear"] = corr_matrix["shear"]["opening"] 

148 corr_matrix.loc["shear", "opening"] = corr_matrix.loc["opening", "shear"] 

149 

150 # stretch 

151 # corr_matrix["stretch"]["stagger"] = stretch.corrwith(stagger, method="pearson") 

152 corr_matrix.loc["stagger", "stretch"] = stretch.corrwith(stagger, method="pearson").values[0] 

153 # corr_matrix["stretch"]["buckle"] = stretch.corrwith(buckle, method=self.circlineal) 

154 corr_matrix.loc["buckle", "stretch"] = stretch.corrwith(buckle, method=self.circlineal).values[0] # type: ignore 

155 # corr_matrix["stretch"]["propel"] = stretch.corrwith(propel, method=self.circlineal) 

156 corr_matrix.loc["propel", "stretch"] = stretch.corrwith(propel, method=self.circlineal).values[0] # type: ignore 

157 # corr_matrix["stretch"]["opening"] = stretch.corrwith(opening, method=self.circlineal) 

158 corr_matrix.loc["opening", "stretch"] = stretch.corrwith(opening, method=self.circlineal).values[0] # type: ignore 

159 # symmetric values 

160 # corr_matrix["stagger"]["stretch"] = corr_matrix["stretch"]["stagger"] 

161 corr_matrix.loc["stretch", "stagger"] = corr_matrix.loc["stagger", "stretch"] 

162 # corr_matrix["buckle"]["stretch"] = corr_matrix["stretch"]["buckle"] 

163 corr_matrix.loc["stretch", "buckle"] = corr_matrix.loc["buckle", "stretch"] 

164 # corr_matrix["propel"]["stretch"] = corr_matrix["stretch"]["propel"] 

165 corr_matrix.loc["stretch", "propel"] = corr_matrix.loc["propel", "stretch"] 

166 # corr_matrix["opening"]["stretch"] = corr_matrix["stretch"]["opening"] 

167 corr_matrix.loc["stretch", "opening"] = corr_matrix.loc["opening", "stretch"] 

168 

169 # stagger 

170 # corr_matrix["stagger"]["buckle"] = stagger.corrwith(buckle, method=self.circlineal) 

171 corr_matrix.loc["buckle", "stagger"] = stagger.corrwith(buckle, method=self.circlineal).values[0] # type: ignore 

172 # corr_matrix["stagger"]["propel"] = stagger.corrwith(propel, method=self.circlineal) 

173 corr_matrix.loc["propel", "stagger"] = stagger.corrwith(propel, method=self.circlineal).values[0] # type: ignore 

174 # corr_matrix["stagger"]["opening"] = stagger.corrwith(opening, method=self.circlineal) 

175 corr_matrix.loc["opening", "stagger"] = stagger.corrwith(opening, method=self.circlineal).values[0] # type: ignore 

176 # symmetric values 

177 # corr_matrix["buckle"]["stagger"] = corr_matrix["stagger"]["buckle"] 

178 corr_matrix.loc["stagger", "buckle"] = corr_matrix.loc["buckle", "stagger"] 

179 # corr_matrix["propel"]["stagger"] = corr_matrix["stagger"]["propel"] 

180 corr_matrix.loc["stagger", "propel"] = corr_matrix.loc["propel", "stagger"] 

181 # corr_matrix["opening"]["stagger"] = corr_matrix["stagger"]["opening"] 

182 corr_matrix.loc["stagger", "opening"] = corr_matrix.loc["opening", "stagger"] 

183 

184 # buckle 

185 # corr_matrix["buckle"]["propel"] = buckle.corrwith(propel, method=self.circular) 

186 corr_matrix.loc["propel", "buckle"] = buckle.corrwith(propel, method=self.circular).values[0] # type: ignore 

187 # corr_matrix["buckle"]["opening"] = buckle.corrwith(opening, method=self.circular) 

188 corr_matrix.loc["opening", "buckle"] = buckle.corrwith(opening, method=self.circular).values[0] # type: ignore 

189 # symmetric values 

190 # corr_matrix["propel"]["buckle"] = corr_matrix["buckle"]["propel"] 

191 corr_matrix.loc["buckle", "propel"] = corr_matrix.loc["propel", "buckle"] 

192 # corr_matrix["opening"]["buckle"] = corr_matrix["buckle"]["opening"] 

193 corr_matrix.loc["buckle", "opening"] = corr_matrix.loc["opening", "buckle"] 

194 

195 # propel 

196 # corr_matrix["propel"]["opening"] = propel.corrwith(opening, method=self.circular) 

197 corr_matrix.loc["opening", "propel"] = propel.corrwith(opening, method=self.circular).values[0] # type: ignore 

198 # symmetric values 

199 # corr_matrix["opening"]["propel"] = corr_matrix["propel"]["opening"] 

200 corr_matrix.loc["propel", "opening"] = corr_matrix.loc["opening", "propel"] 

201 

202 # save csv data 

203 corr_matrix.to_csv(self.stage_io_dict["out"]["output_csv_path"]) 

204 

205 # create heatmap 

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

207 axs.pcolor(corr_matrix) 

208 # Loop over data dimensions and create text annotations. 

209 for i in range(len(corr_matrix)): 

210 for j in range(len(corr_matrix)): 

211 axs.text( 

212 j+.5, 

213 i+.5, 

214 f"{corr_matrix[coordinates[j]].loc[coordinates[i]]:.2f}", 

215 ha="center", 

216 va="center", 

217 color="w") 

218 axs.set_xticks([i + 0.5 for i in range(len(corr_matrix))]) 

219 axs.set_xticklabels(corr_matrix.columns, rotation=90) 

220 axs.set_yticks([i+0.5 for i in range(len(corr_matrix))]) 

221 axs.set_yticklabels(corr_matrix.index) 

222 axs.set_title( 

223 "Helical Parameter Correlation " 

224 f"for Base Pair Step \'{self.base}\'") 

225 fig.tight_layout() 

226 fig.savefig( 

227 self.stage_io_dict['out']['output_jpg_path'], 

228 format="jpg") 

229 plt.close() 

230 

231 # Copy files to host 

232 self.copy_to_host() 

233 

234 # Remove temporary file(s) 

235 self.remove_tmp_files() 

236 

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

238 

239 return 0 

240 

241 def get_corr_method(self, corrtype1, corrtype2): 

242 if corrtype1 == "circular" and corrtype2 == "linear": 

243 method = self.circlineal 

244 if corrtype1 == "linear" and corrtype2 == "circular": 

245 method = self.circlineal 

246 elif corrtype1 == "circular" and corrtype2 == "circular": 

247 method = self.circular 

248 else: 

249 method = "pearson" 

250 return method 

251 

252 @staticmethod 

253 def circular(x1, x2): 

254 x1 = x1 * np.pi / 180 

255 x2 = x2 * np.pi / 180 

256 diff_1 = np.sin(x1 - x1.mean()) 

257 diff_2 = np.sin(x2 - x2.mean()) 

258 num = (diff_1 * diff_2).sum() 

259 den = np.sqrt((diff_1 ** 2).sum() * (diff_2 ** 2).sum()) 

260 return num / den 

261 

262 @staticmethod 

263 def circlineal(x1, x2): 

264 x2 = x2 * np.pi / 180 

265 rc = np.corrcoef(x1, np.cos(x2))[1, 0] 

266 rs = np.corrcoef(x1, np.sin(x2))[1, 0] 

267 rcs = np.corrcoef(np.sin(x2), np.cos(x2))[1, 0] 

268 num = (rc ** 2) + (rs ** 2) - 2 * rc * rs * rcs 

269 den = 1 - (rcs ** 2) 

270 correlation = np.sqrt(num / den) 

271 if np.corrcoef(x1, x2)[1, 0] < 0: 

272 correlation *= -1 

273 return correlation 

274 

275 

276def intrahpcorr( 

277 input_filename_shear: str, input_filename_stretch: str, 

278 input_filename_stagger: str, input_filename_buckle: str, 

279 input_filename_propel: str, input_filename_opening: str, 

280 output_csv_path: str, output_jpg_path: str, 

281 properties: Optional[dict] = None, **kwargs) -> int: 

282 """Create :class:`IntraHelParCorrelation <intrabp_correlations.intrahpcorr.IntraHelParCorrelation>` class and 

283 execute the :meth:`launch() <intrabp_correlations.intrahpcorr.IntraHelParCorrelation.launch>` method.""" 

284 return IntraHelParCorrelation(**dict(locals())).launch() 

285 

286 

287intrahpcorr.__doc__ = IntraHelParCorrelation.__doc__ 

288main = IntraHelParCorrelation.get_main(intrahpcorr, "Load helical parameter file and save base data individually.") 

289 

290if __name__ == '__main__': 

291 main()