Coverage for bbconf/bbagent.py: 72%

1import datetime

2import logging

3import statistics

4import threading

5from functools import cached_property

6from pathlib import Path

8import numpy as np

9from cachetools import TTLCache

10from sqlalchemy.engine import ScalarResult

11from sqlalchemy.orm import Session

12from sqlalchemy.sql import and_, distinct, func, or_, select

14from bbconf.config_parser.bedbaseconfig import BedBaseConfig

15from bbconf.db_utils import (

16 Bed,

17 BedMetadata,

18 BedSets,

19 BedStats,

20 Files,

21 GeoGsmStatus,

22 License,

23 ReferenceGenome,

24 UsageBedMeta,

25 UsageBedSetMeta,

26 UsageFiles,

27 UsageSearch,

28)

29from bbconf.models.base_models import (

30 AllFilesInfo,

31 BinValues,

32 FileInfo,

33 FileStats,

34 GEOStatistics,

35 StatsReturn,

36 UsageItem,

37 UsageModel,

38 UsageResponse,

39 UsageStats,

40)

41from bbconf.modules.analysis_files import BedAgentAnalysisFile

42from bbconf.modules.bedfiles import BedAgentBedFile

43from bbconf.modules.bedsets import BedAgentBedSet

44from bbconf.modules.objects import BBObjects

45from bbconf.modules.snapshots import BedAgentSnapshot

47from .const import PKG_NAME

49_LOGGER = logging.getLogger(PKG_NAME)

52class BedBaseAgent:

53 def __init__(

54 self,

55 config: Path | str,

56 init_ml: bool = True,

57 ):

58 """

59 Initialize connection to the pep_db database. You can use the basic connection parameters

60 or libpq connection string.

62 Args:

63 config: Path to the configuration file.

64 init_ml: Initialize ML models for search (default: True).

65 """

67 self.config = BedBaseConfig(config, init_ml)

69 self._bed = BedAgentBedFile(self.config, self)

70 self._bedset = BedAgentBedSet(self.config)

71 self._objects = BBObjects(self.config)

72 self._snapshot = BedAgentSnapshot(self.config)

73 self._analysis_files = BedAgentAnalysisFile(self.config)

75 # get_stats() runs three uncached COUNT queries on the multi-hundred-

76 # thousand-row bed table and is called on hot paths (the stats endpoint

77 # plus the neighbours/list/search result builders). Cache the result

78 # with a TTL so those paths do not hit the database on every request.

79 # The lock guards the cache dict only, never the DB query itself.

80 self._stats_cache = TTLCache(maxsize=1, ttl=3600)

81 self._stats_lock = threading.Lock()

83 @property

84 def bed(self) -> BedAgentBedFile:

85 return self._bed

87 @property

88 def bedset(self) -> BedAgentBedSet:

89 return self._bedset

91 @property

92 def objects(self) -> BBObjects:

93 return self._objects

95 @property

96 def snapshot(self) -> BedAgentSnapshot:

97 return self._snapshot

99 @property

100 def analysis_files(self) -> BedAgentAnalysisFile:

101 return self._analysis_files

102

103 def __repr__(self) -> str:

104 repr = f"BedBaseAgent(config={self.config})"

105 repr += f"\n{self.bed}"

106 repr += f"\n{self.bedset}"

107 repr += f"\n{self.objects}"

108 return repr

109

110 def get_stats(self) -> StatsReturn:

111 """

112 Get statistics for a bed file.

113

114 The result is cached with a TTL because this runs three COUNT queries

115 against the large bed table and is called on hot API paths.

116

117 Returns:

118 Statistics.

119 """

120 with self._stats_lock:

121 cached = self._stats_cache.get("stats")

122 if cached is not None:

123 return cached

124

125 with Session(self.config.db_engine.engine) as session:

126 number_of_bed = session.execute(select(func.count(Bed.id))).one()[0]

127 number_of_bedset = session.execute(select(func.count(BedSets.id))).one()[0]

128

129 number_of_genomes = session.execute(

130 select(func.count(distinct(Bed.genome_alias)))

131 ).one()[0]

132

133 stats = StatsReturn(

134 bedfiles_number=number_of_bed,

135 bedsets_number=number_of_bedset,

136 genomes_number=number_of_genomes,

137 )

138

139 with self._stats_lock:

140 self._stats_cache["stats"] = stats

141

142 return stats

143

144 def get_detailed_stats(self, concise: bool = False) -> FileStats:

145 """

146 Get comprehensive statistics for all bed files.

147

148 Args:

149 concise: If True, return only top 20 items for each category.

150

151 Returns:

152 FileStats object containing detailed statistics.

153 """

154

155 _LOGGER.info("Getting detailed statistics for all bed files")

156

157 numeric_stats_statement = (

158 select(

159 BedStats.number_of_regions,

160 BedStats.mean_region_width,

161 Files.size,

162 )

163 .select_from(Bed)

164 .join(BedStats, BedStats.id == Bed.id)

165 .join(Files, Files.bedfile_id == Bed.id)

166 .where(

167 Files.name == "bed_file",

168 BedStats.number_of_regions.is_not(None),

169 BedStats.mean_region_width.is_not(None),

170 Files.size.is_not(None),

171 )

172 )

173

174 with Session(self.config.db_engine.engine) as session:

175 bed_compliance = {

176 f[0]: f[1]

177 for f in session.execute(

178 select(Bed.bed_compliance, func.count(Bed.bed_compliance))

179 .where(Bed.bed_compliance.is_not(None))

180 .group_by(Bed.bed_compliance)

181 .order_by(func.count(Bed.bed_compliance).desc())

182 ).all()

183 }

184 data_format = {

185 f[0]: f[1]

186 for f in session.execute(

187 select(Bed.data_format, func.count(Bed.data_format))

188 .where(Bed.data_format.is_not(None))

189 .group_by(Bed.data_format)

190 .order_by(func.count(Bed.data_format).desc())

191 ).all()

192 }

193 file_genomes = {

194 f[0]: f[1]

195 for f in session.execute(

196 select(Bed.genome_alias, func.count(Bed.genome_alias))

197 .where(Bed.genome_alias.is_not(None))

198 .group_by(Bed.genome_alias)

199 .order_by(func.count(Bed.genome_alias).desc())

200 ).all()

201 }

202 file_organism = {

203 f[0]: f[1]

204 for f in session.execute(

205 select(

206 BedMetadata.species_name, func.count(BedMetadata.species_name)

207 )

208 .where(BedMetadata.species_name.is_not(None))

209 .group_by(BedMetadata.species_name)

210 .order_by(func.count(BedMetadata.species_name).desc())

211 ).all()

212 }

213 file_assay = {

214 f[0]: f[1]

215 for f in session.execute(

216 select(BedMetadata.assay, func.count(BedMetadata.assay))

217 .where(BedMetadata.assay.is_not(None))

218 .group_by(BedMetadata.assay)

219 .order_by(func.count(BedMetadata.assay).desc())

220 ).all()

221 }

222 cell_line = {

223 f[0]: f[1]

224 for f in session.execute(

225 select(BedMetadata.cell_line, func.count(BedMetadata.cell_line))

226 .where(BedMetadata.cell_line.is_not(None))

227 .group_by(BedMetadata.cell_line)

228 .order_by(func.count(BedMetadata.cell_line).desc())

229 ).all()

230 }

231

232 bed_comments = self._stats_comments(session)

233 geo_status = self._stats_geo_status(session)

234

235 numeric_rows = session.execute(numeric_stats_statement).all()

236

237 geo_stats = self._get_geo_stats(session)

238

239 slice_value = 20

240

241 number_of_regions = [row[0] for row in numeric_rows]

242 list_mean_width = [row[1] for row in numeric_rows]

243 list_file_size = [row[2] for row in numeric_rows]

244

245 number_of_regions_bins = self._bin_number_of_regions(number_of_regions)

246 list_mean_width_bins = self._bin_mean_region_width(list_mean_width)

247 list_file_size_bins = self._bin_file_size(list_file_size)

248

249 if concise:

250 bed_compliance_concise = dict(list(bed_compliance.items())[0:slice_value])

251 bed_compliance_concise["other"] = sum(

252 list(bed_compliance.values())[slice_value:]

253 )

254 if "" in bed_compliance_concise:

255 bed_compliance_concise["other"] = (

256 bed_compliance_concise["other"] + bed_compliance_concise[""]

257 )

258 bed_compliance_concise.pop("")

259

260 file_genomes_concise = dict(list(file_genomes.items())[0:slice_value])

261 file_genomes_concise["other"] = sum(

262 list(file_genomes.values())[slice_value:]

263 ) + file_genomes.get("other", 0)

264 if "" in file_genomes_concise:

265 file_genomes_concise["other"] = (

266 file_genomes_concise["other"] + file_genomes_concise[""]

267 )

268 file_genomes_concise.pop("")

269

270 file_organism_concise = dict(list(file_organism.items())[0:slice_value])

271 file_organism_concise["other"] = sum(

272 list(file_organism.values())[slice_value:]

273 ) + file_organism.get("other", 0)

274 if "" in file_organism_concise:

275 file_organism_concise["other"] = (

276 file_organism_concise["other"] + file_organism_concise[""]

277 )

278 file_organism_concise.pop("")

279 file_assay_concise = dict(list(file_assay.items())[0:slice_value])

280 file_assay_concise["other"] = sum(

281 list(file_assay.values())[slice_value:]

282 ) + file_assay.get("other", 0)

283 if "" in file_assay_concise:

284 file_assay_concise["other"] = (

285 file_assay_concise["other"] + file_assay_concise[""]

286 )

287 file_assay_concise.pop("")

288 if "OTHER" in file_assay_concise:

289 file_assay_concise["other"] = (

290 file_assay_concise["other"] + file_assay_concise["OTHER"]

291 )

292 file_assay_concise.pop("OTHER")

293

294 cell_line_concise = dict(list(cell_line.items())[0:slice_value])

295 cell_line_concise["other"] = sum(

296 list(cell_line.values())[slice_value:]

297 ) + cell_line.get("other", 0)

298 if "" in cell_line_concise:

299 cell_line_concise["other"] = (

300 cell_line_concise["other"] + cell_line_concise[""]

301 )

302 cell_line_concise.pop("")

303

304 return FileStats(

305 data_format=data_format,

306 bed_compliance=bed_compliance_concise,

307 file_genome=file_genomes_concise,

308 file_organism=file_organism_concise,

309 file_assay=file_assay_concise,

310 cell_line=cell_line_concise,

311 bed_comments=bed_comments,

312 geo_status=geo_status,

313 mean_region_width=list_mean_width_bins,

314 file_size=list_file_size_bins,

315 number_of_regions=number_of_regions_bins,

316 geo=geo_stats,

317 )

318

319 return FileStats(

320 data_format=data_format,

321 bed_compliance=bed_compliance,

322 file_genome=file_genomes,

323 file_organism=file_organism,

324 file_assay=file_assay,

325 cell_line=cell_line,

326 bed_comments=bed_comments,

327 geo_status=geo_status,

328 mean_region_width=list_mean_width_bins,

329 file_size=list_file_size_bins,

330 number_of_regions=number_of_regions_bins,

331 geo=geo_stats,

332 )

333

334 def get_detailed_usage(self) -> UsageStats:

335 """

336 Get detailed usage statistics for the bedbase platform.

337 This method will only return top 20 items for each category.

338

339 Returns:

340 UsageStats object containing detailed usage statistics.

341 """

342

343 _LOGGER.info("Getting detailed usage statistics.")

344

345 with Session(self.config.db_engine.engine) as session:

346 bed_metadata = {

347 f[0]: f[1]

348 for f in session.execute(

349 select(UsageBedMeta.bed_id, func.sum(UsageBedMeta.count))

350 .group_by(UsageBedMeta.bed_id)

351 .order_by(func.sum(UsageBedMeta.count).desc())

352 .limit(20)

353 ).all()

354 }

355 bedset_metadata = {

356 f[0]: f[1]

357 for f in session.execute(

358 select(UsageBedSetMeta.bedset_id, func.sum(UsageBedSetMeta.count))

359 .group_by(UsageBedSetMeta.bedset_id)

360 .order_by(func.sum(UsageBedSetMeta.count).desc())

361 .limit(20)

362 ).all()

363 }

364

365 bed_search_terms = {

366 f[0]: f[1]

367 for f in session.execute(

368 select(UsageSearch.query, func.sum(UsageSearch.count))

369 .where(UsageSearch.type == "bed")

370 .group_by(UsageSearch.query)

371 .order_by(func.sum(UsageSearch.count).desc())

372 .limit(20)

373 ).all()

374 if f[0]

375 }

376

377 bedset_search_terms = {

378 f[0]: f[1]

379 for f in session.execute(

380 select(UsageSearch.query, func.sum(UsageSearch.count))

381 .where(UsageSearch.type == "bedset")

382 .group_by(UsageSearch.query)

383 .order_by(func.sum(UsageSearch.count).desc())

384 .limit(20)

385 ).all()

386 if f[0]

387 }

388

389 bed_downloads = {

390 f[0].split("/")[-1].split(".")[0]: f[1]

391 for f in session.execute(

392 select(UsageFiles.file_path, func.sum(UsageFiles.count))

393 .group_by(UsageFiles.file_path)

394 .order_by(func.sum(UsageFiles.count).desc())

395 .limit(20)

396 ).all()

397 }

398

399 return UsageStats(

400 bed_metadata=bed_metadata,

401 bedset_metadata=bedset_metadata,

402 bed_search_terms=bed_search_terms,

403 bedset_search_terms=bedset_search_terms,

404 bed_downloads=bed_downloads,

405 )

406

407 def get_usage(

408 self,

409 event_type: str,

410 date_from: datetime.datetime | None = None,

411 date_to: datetime.datetime | None = None,

412 limit: int = 100,

413 offset: int = 0,

414 ) -> UsageResponse:

415 """

416 Get per-key usage counts for a single event type, optionally restricted

417 to a date range.

418

419 Unlike get_detailed_usage (which sums across all time and is hardcoded to

420 the top 20 keys), this method groups by key and sums the counts of the

421 time-bucketed rows whose window overlaps the requested date range, and

422 supports pagination. It is read-only and does not modify any usage data.

423

424 Args:

425 event_type: One of "files", "bed_meta", "bedset_meta", "bed_search",

426 "bedset_search".

427 date_from: Include only rows whose window ends on or after this date

428 (date_to >= date_from). None means no lower bound.

429 date_to: Include only rows whose window starts on or before this date

430 (date_from <= date_to). None means no upper bound.

431 limit: Maximum number of keys to return (>= 1).

432 offset: Number of keys to skip (for pagination).

433

434 Returns:

435 UsageResponse with per-key counts ordered by count descending.

436 """

437 # model, key column, and (for UsageSearch) the `type` filter value.

438 usage_map = {

439 "files": (UsageFiles, UsageFiles.file_path, None),

440 "bed_meta": (UsageBedMeta, UsageBedMeta.bed_id, None),

441 "bedset_meta": (UsageBedSetMeta, UsageBedSetMeta.bedset_id, None),

442 "bed_search": (UsageSearch, UsageSearch.query, "bed"),

443 "bedset_search": (UsageSearch, UsageSearch.query, "bedset"),

444 }

445

446 if event_type not in usage_map:

447 raise ValueError(

448 f"Invalid event_type '{event_type}'. Must be one of: "

449 f"{', '.join(usage_map)}."

450 )

451 if limit < 1 or offset < 0:

452 raise ValueError("limit must be >= 1 and offset >= 0.")

453

454 model, key_column, search_type = usage_map[event_type]

455

456 _LOGGER.info(f"Getting usage statistics for event type '{event_type}'.")

457

458 statement = (

459 select(key_column, func.sum(model.count))

460 .group_by(key_column)

461 .order_by(func.sum(model.count).desc())

462 )

463

464 statement = statement.where(key_column.isnot(None), key_column != "")

465 if search_type is not None:

466 statement = statement.where(model.type == search_type)

467 if date_from is not None:

468 statement = statement.where(model.date_to >= date_from)

469 if date_to is not None:

470 statement = statement.where(model.date_from <= date_to)

471

472 statement = statement.limit(limit).offset(offset)

473

474 # bed downloads are keyed by file path; shorten to the bed id, matching

475 # get_detailed_usage.

476 shorten_key = event_type == "files"

477

478 with Session(self.config.db_engine.engine) as session:

479 results = [

480 UsageItem(

481 key=(f[0].split("/")[-1].split(".")[0] if shorten_key else f[0]),

482 count=f[1],

483 )

484 for f in session.execute(statement).all()

485 ]

486

487 return UsageResponse(

488 event_type=event_type,

489 date_from=date_from,

490 date_to=date_to,

491 limit=limit,

492 offset=offset,

493 count=len(results),

494 results=results,

495 )

496

497 def get_list_genomes(self) -> list[str]:

498 """

499 Get list of genomes from the database.

500

501 Returns:

502 List of genomes.

503 """

504 statement = (

505 select(Bed.genome_alias)

506 .group_by(Bed.genome_alias)

507 .order_by(func.count(Bed.genome_alias).desc())

508 )

509 with Session(self.config.db_engine.engine) as session:

510 genomes = session.execute(statement).all()

511 return [result[0] for result in genomes if result[0]]

512

513 def get_list_assays(self) -> list[str]:

514 """

515 Get list of assays from the database.

516

517 Returns:

518 List of assays.

519 """

520

521 with Session(self.config.db_engine.engine) as session:

522 statement = (

523 select(BedMetadata.assay)

524 .group_by(BedMetadata.assay)

525 .order_by(func.count(BedMetadata.assay).desc())

526 )

527 results = session.execute(statement).all()

528 return [result[0] for result in results if result[0]]

529

530 @cached_property

531 def list_of_licenses(self) -> list[str]:

532 """

533 Get list of licenses from the database.

534

535 Returns:

536 List of licenses.

537 """

538 statement = select(License.id)

539 with Session(self.config.db_engine.engine) as session:

540 licenses = session.execute(statement).all()

541 return [result[0] for result in licenses]

542

543 def add_usage(self, stats: UsageModel) -> None:

544

545 with Session(self.config.db_engine.engine) as session:

546 # FILES USAGE

547 reported_items_files = session.scalars(

548 select(UsageFiles).where(UsageFiles.date_to > func.now())

549 )

550 reported_dict_files = {

551 item.file_path: item for item in reported_items_files

552 }

553

554 for key, value in stats.files.items():

555 if key in reported_dict_files:

556 reported_dict_files[key].count += value

557 else:

558 new_stats = UsageFiles(

559 file_path=key,

560 count=value,

561 date_from=stats.date_from,

562 date_to=stats.date_to,

563 )

564 session.add(new_stats)

565

566 # METADATA USAGE

567 reported_items_metadata = session.scalars(

568 select(UsageBedMeta).where(UsageBedMeta.date_to > func.now())

569 )

570

571 reported_dict_metadata = {

572 item.bed_id: item for item in reported_items_metadata

573 }

574

575 for key, value in stats.bed_meta.items():

576 if key in reported_dict_metadata:

577 reported_dict_metadata[key].count += value

578 else:

579 new_stats = UsageBedMeta(

580 bed_id=key,

581 count=value,

582 date_from=stats.date_from,

583 date_to=stats.date_to,

584 )

585 session.add(new_stats)

586

587 # BEDSET METADATA USAGE

588 reported_items_bedset_metadata = session.scalars(

589 select(UsageBedSetMeta).where(UsageBedSetMeta.date_to > func.now())

590 )

591

592 reported_dict_bedset_metadata = {

593 item.bedset_id: item for item in reported_items_bedset_metadata

594 }

595

596 for key, value in stats.bedset_meta.items():

597 if key in reported_dict_bedset_metadata:

598 reported_dict_bedset_metadata[key].count += value

599 else:

600 new_stats = UsageBedSetMeta(

601 bedset_id=key,

602 count=value,

603 date_from=stats.date_from,

604 date_to=stats.date_to,

605 )

606 session.add(new_stats)

607

608 # SEARCH USAGE

609

610 reported_items_bed_search = session.scalars(

611 select(UsageSearch).where(

612 UsageSearch.type == "bed", UsageSearch.date_to > func.now()

613 )

614 )

615 reported_dict_bed_search = {

616 item.query: item for item in reported_items_bed_search

617 }

618

619 for key, value in stats.bed_search.items():

620 if key in reported_dict_bed_search:

621 reported_dict_bed_search[key].count += value

622 else:

623 new_stats = UsageSearch(

624 query=key,

625 count=value,

626 type="bed",

627 date_from=stats.date_from,

628 date_to=stats.date_to,

629 )

630 session.add(new_stats)

631

632 # SEARCH BEDSET USAGE

633 reporeted_items_bedset_search = session.scalars(

634 select(UsageSearch).where(

635 UsageSearch.type == "bedset", UsageSearch.date_to > func.now()

636 )

637 )

638 reported_dict_bedset_search = {

639 item.query: item for item in reporeted_items_bedset_search

640 }

641

642 for key, value in stats.bedset_search.items():

643 if key in reported_dict_bedset_search:

644 reported_dict_bedset_search[key].count += value

645 else:

646 new_stats = UsageSearch(

647 query=key,

648 count=value,

649 type="bedset",

650 date_from=stats.date_from,

651 date_to=stats.date_to,

652 )

653 session.add(new_stats)

654

655 session.commit()

656

657 def _stats_comments(self, sa_session: Session) -> dict[str, int]:

658 """

659 Get statistics about comments that are present in bed files.

660

661 Args:

662 sa_session: SQLAlchemy session.

663

664 Returns:

665 Dict mapping comment type to count.

666 """

667 _LOGGER.info("Analyzing bed table for comments in bed files...")

668

669 total_statement = select(func.count(Bed.id)).where(Bed.header.is_not(None))

670 correct_statement = select(func.count(Bed.id)).where(Bed.header.like("#%"))

671 track_name_statement = select(func.count(Bed.id)).where(

672 Bed.header.like("track name%")

673 )

674 track_type_statement = select(func.count(Bed.id)).where(

675 Bed.header.like("track type%")

676 )

677 browser_statement = select(func.count(Bed.id)).where(

678 Bed.header.like("browser%")

679 )

680

681 total_bed_comments = (

682 sa_session.execute(total_statement).one_or_none() or (0,)

683 )[0]

684 correct_bed_comments = (

685 sa_session.execute(correct_statement).one_or_none() or (0,)

686 )[0]

687 track_name_comments = (

688 sa_session.execute(track_name_statement).one_or_none() or (0,)

689 )[0]

690 track_type_comments = (

691 sa_session.execute(track_type_statement).one_or_none() or (0,)

692 )[0]

693 browser_comments = (

694 sa_session.execute(browser_statement).one_or_none() or (0,)

695 )[0]

696

697 header_comments = (

698 total_bed_comments

699 - correct_bed_comments

700 - track_name_comments

701 - track_type_comments

702 - browser_comments

703 )

704

705 return {

706 "correct_bed_comments": correct_bed_comments,

707 "track_name_comments": track_name_comments,

708 "track_type_comments": track_type_comments,

709 "browser_comments": browser_comments,

710 "header_comments": header_comments,

711 }

712

713 def _stats_geo_status(self, sa_session: Session) -> dict[str, int]:

714 """

715 Get statistics about status of GEO bed file processing.

716

717 Args:

718 sa_session: SQLAlchemy session.

719

720 Returns:

721 Dict mapping status type to count.

722 """

723

724 success_statement = select(

725 func.count(distinct(GeoGsmStatus.sample_name))

726 ).where(GeoGsmStatus.status == "SUCCESS")

727

728 failed_count_statement = select(

729 func.count(distinct(GeoGsmStatus.sample_name))

730 ).where(and_(GeoGsmStatus.status == "FAIL"))

731 mean_region_width_statement = select(

732 func.count(distinct(GeoGsmStatus.sample_name))

733 ).where(

734 or_(

735 GeoGsmStatus.error.like("%Initial QC failed for%"),

736 GeoGsmStatus.error.like("%Quality control failed%"),

737 )

738 )

739 size_greater_then_statement = select(

740 func.count(distinct(GeoGsmStatus.sample_name))

741 ).where(GeoGsmStatus.error.like("%File size is too big.%"))

742

743 success_count = (sa_session.execute(success_statement).one_or_none() or (0,))[0]

744 failed_count = (

745 sa_session.execute(failed_count_statement).one_or_none() or (0,)

746 )[0]

747

748 mean_region_width_count = (

749 sa_session.execute(mean_region_width_statement).one_or_none() or (0,)

750 )[0]

751 size_greater_than_count = (

752 sa_session.execute(size_greater_then_statement).one_or_none() or (0,)

753 )[0]

754 corrupted_file_count = (

755 failed_count - mean_region_width_count - size_greater_than_count

756 )

757

758 return {

759 "success": success_count,

760 "mean_region_width_lower_10": mean_region_width_count,

761 "size_greater_20": size_greater_than_count,

762 "corrupted_files": corrupted_file_count,

763 }

764

765 def bed_files_info(self) -> AllFilesInfo:

766 """

767 Get information about all bed files in bedbase.

768

769 Returns:

770 AllFilesInfo containing total count and list of file info objects.

771 """

772

773 all_files_statement = (

774 select(

775 Bed.id,

776 Bed.bed_compliance,

777 Bed.data_format,

778 BedStats.mean_region_width,

779 Files.size,

780 BedStats.number_of_regions,

781 )

782 .join(BedStats, BedStats.id == Bed.id) # Explicit join condition

783 .join(Files, Files.bedfile_id == Bed.id) # Explicit join condition

784 .where(Files.name == "bed_file")

785 )

786

787 results = []

788 error_list = []

789 with Session(self.config.db_engine.engine) as session:

790 sql_result = session.execute(all_files_statement).all()

791

792 for single_bed in sql_result:

793 try:

794 results.append(

795 FileInfo(

796 id=single_bed[0],

797 bed_compliance=single_bed[1],

798 data_format=single_bed[2],

799 mean_region_width=single_bed[3],

800 file_size=single_bed[4],

801 number_of_regions=single_bed[5],

802 )

803 )

804 except Exception:

805 error_list.append(single_bed[0])

806

807 _LOGGER.info(

808 f"Number of bed records with unknown regions and region width {len(error_list)}"

809 )

810 return AllFilesInfo(

811 total=len(results),

812 files=results,

813 )

814

815 def _bin_number_of_regions(self, number_of_regions: list) -> BinValues:

816 """

817 Create bins for number of regions in bed files.

818

819 Args:

820 number_of_regions: List of number of regions in bed files.

821

822 Returns:

823 BinValues object containing bins and values.

824 """

825

826 max_value_threshold = 400_000 # set a threshold for maximum value to avoid outliers in the histogram

827

828 filtered_number_of_regions = [

829 x if x <= max_value_threshold else max_value_threshold + 1

830 for x in number_of_regions

831 ]

832

833 n_region_counts, n_region_bin_edges = np.histogram(

834 filtered_number_of_regions, bins=100

835 )

836 n_region_counts = n_region_counts.astype(int).tolist()

837 n_region_bin_edges = n_region_bin_edges.astype(int).tolist()

838

839 return BinValues(

840 bins=n_region_bin_edges,

841 counts=n_region_counts,

842 mean=round(statistics.mean(number_of_regions), 2)

843 if number_of_regions

844 else 0,

845 median=(

846 round(statistics.median(number_of_regions), 2)

847 if number_of_regions

848 else 0

849 ),

850 )

851

852 def _bin_mean_region_width(self, mean_region_widths: list) -> BinValues:

853 """

854 Create bins for number of regions in bed files.

855

856 Args:

857 mean_region_widths: List of mean region widths in bed files.

858

859 Returns:

860 BinValues object containing bins and values.

861 """

862

863 max_value_threshold = 5_000 # set a threshold for maximum value to avoid outliers in the histogram

864

865 filtered_mean_region_widths = [

866 x if x <= max_value_threshold else max_value_threshold + 1

867 for x in mean_region_widths

868 ]

869

870 mean_reg_width_counts, mean_reg_width_bin_edges = np.histogram(

871 filtered_mean_region_widths, bins=100

872 )

873 mean_reg_width_counts = mean_reg_width_counts.tolist()

874 mean_reg_width_bin_edges = mean_reg_width_bin_edges.tolist()

875

876 return BinValues(

877 bins=mean_reg_width_bin_edges,

878 counts=mean_reg_width_counts,

879 mean=(

880 round(statistics.mean(mean_region_widths), 2)

881 if mean_region_widths

882 else 0

883 ),

884 median=(

885 round(statistics.median(mean_region_widths), 2)

886 if mean_region_widths

887 else 0

888 ),

889 )

890

891 def _bin_file_size(self, list_file_size: list) -> BinValues:

892 """

893 Create bins for number of regions in bed files.

894

895 Args:

896 list_file_size: List of bed file sizes in bytes.

897

898 Returns:

899 BinValues object containing bins and values.

900 """

901

902 max_value_threshold = 10 * 1024 * 1024

903

904 filtered_list_file_size = [

905 x for x in list_file_size if x <= max_value_threshold

906 ]

907

908 filtered_list_file_size = [x / (1024 * 1024) for x in filtered_list_file_size]

909

910 file_size_counts, file_size_bin_edges = np.histogram(

911 filtered_list_file_size, bins=100

912 )

913 file_size_counts = file_size_counts.astype(int).tolist()

914 file_size_bin_edges = file_size_bin_edges.astype(float).tolist()

915

916 return BinValues(

917 bins=file_size_bin_edges,

918 counts=file_size_counts,

919 mean=(

920 round(statistics.mean(filtered_list_file_size), 2)

921 if filtered_list_file_size

922 else 0

923 ),

924 median=(

925 round(statistics.median(filtered_list_file_size), 2)

926 if filtered_list_file_size

927 else 0

928 ),

929 )

930

931 def _get_geo_stats(self, sa_session: Session) -> GEOStatistics:

932 """

933 Get GEO statistics for the bedbase platform.

934

935 Returns:

936 GEOStatistics.

937 """

938

939 _LOGGER.info("Getting GEO statistics.")

940

941 statement = select(

942 GeoGsmStatus.bed_id,

943 GeoGsmStatus.source_submission_date,

944 GeoGsmStatus.file_size,

945 ).distinct(GeoGsmStatus.sample_name)

946

947 results = sa_session.execute(statement).all()

948

949 years = []

950 file_sizes = []

951

952 for result in results:

953 if result[1]:

954 years.append(result[1].year)

955 if result[2]:

956 file_sizes.append(result[2])

957

958 years_array = np.array([y for y in years if y is not None])

959 unique_years, counts = np.unique(years_array, return_counts=True)

960

961 cumulative_counts = np.cumsum(counts)

962 unique_years = unique_years.astype(str).tolist()

963

964 years_cumulative = dict(zip(unique_years, cumulative_counts))

965 years_numbers = dict(zip(unique_years, counts))

966

967 MAX_FILE_SIZE = 50 * 1024 * 1024 # 50 MB

968 file_size_filter = [

969 x if x <= MAX_FILE_SIZE else MAX_FILE_SIZE + 1 for x in file_sizes

970 ]

971

972 file_size_filter = [x / (1024 * 1024) for x in file_size_filter]

973 file_size_counts, file_size_bin_edges = np.histogram(file_size_filter, bins=100)

974

975 return GEOStatistics(

976 number_of_files=years_numbers,

977 cumulative_number_of_files=years_cumulative,

978 file_sizes=BinValues(

979 bins=list(file_size_bin_edges),

980 counts=file_size_counts.astype(int).tolist(),

981 mean=round(statistics.mean(file_sizes), 2) if file_sizes else 0,

982 median=round(statistics.median(file_sizes), 2) if file_sizes else 0,

983 ),

984 )

985

986 def get_reference_genomes(self) -> dict[str, str]:

987 """

988 Get mapping of genome aliases to reference genome names.

989

990 Returns:

991 Dict mapping genome_alias to reference_genome_name.

992 """

993

994 genomes = {}

995

996 with Session(self.config.db_engine.engine) as session:

997 results: ScalarResult[ReferenceGenome] = session.scalars(

998 select(ReferenceGenome)

999 )

1000

1001 for genome in results:

1002 genomes[genome.digest] = genome.alias

1003

1004 return genomes