Polypharmacy in individuals with multiple long-term conditions (MLTCs) increases the risk of harmful medication interactions. Identifying specific combinations causing harm is essential for patient safety. This study aims to identify clusters of medications associated with elevated risk of medication-related harm. We linked general practice and hospital data in Wales, focusing on individuals with MLTCs in 2019. Bayesian Profile Regression (BPR) was applied to cluster medicines, conditioned on medication-related harm. Harm was defined using ICD-10 codes for adverse drug events in hospital records within 30 days of a new prescription. Models were adjusted for age, sex, and deprivation. Individual medication profiles utilised Anatomical Therapeutic Chemical (ATC) classifications, with models iteratively refined at more granular ATC levels, to isolate high-risk combinations. Analysis at ATC level 2 - pharmacological/therapeutic subgroup - revealed distinct risk profiles. While the majority of people belonged to low-risk clusters characterised by few prescriptions, two high-risk clusters emerged (harm probabilities: 0.18, 0.19). The first had low probabilities for most medication classes but was defined by distinct combinations of drugs for acid-related disorders, constipation, analgesics, and psycholeptics. Conversely, the highest-risk cluster mirrored lower-risk groups regarding high usage of common medication classes, but was distinguished by specific harmful combinations, such as diuretics with antibacterials. BPR effectively uncovers clinically meaningful medication groupings associated with harm. Linked, population-scale data enables comprehensive analysis of real-world prescribing, supporting targeted prescribing for individuals with MLTCs and polypharmacy. Ongoing work will refine these models using granular ATC classifications to pinpoint specific high-risk medication combinations.