Subsurface thermal gradients along the Queen Charlotte plate boundary from a regional mapping of bottom simulating reflectors

The thermal structure of the lithosphere affects plate boundary evolution and seismogenesis, but obtaining extensive constraints on the subsurface thermal gradient can be challenging. Bottom simulating reflectors (BSRs), the expression of the base of the gas hydrate stability zone (GHSZ) in seismic reflection data, can be used to calculate subsurface thermal gradients in the depth range of the GHSZ (typically within ∼1 km of the seafloor). The Queen Charlotte plate boundary (QCPB), a transpressive system separating the North American and Pacific plates, has hosted multiple magnitude 7+ earthquakes and is undergoing subduction initiation in its southern segment. To investigate the thermal state of the QCPB, we map BSRs within a 2021 regional multichannel seismic reflection dataset and use them to calculate subsurface thermal gradients. BSRs are abundant in the southern half of the QCPB, where compressional deformation dominates, and sparse in the north, found only overlying isolated thrust faults. We find an along-strike thermal gradient trend that corresponds to the increasing age of the Pacific plate. Deviations from the along-strike trend are likely related to focused fluid flow, topography, and erosion. Perpendicular to the QCPB we find a rapid decrease in thermal gradient landward of the deformation front similar to subduction zone prisms and a thermal gradient increase approaching the strike-slip Queen Charlotte Fault, which may reflect fluid flow or shear heating. Differences between previous thermal models and the BSR-derived constraints suggest that the temperature of the Haida Gwaii thrust, and therefore the size of the seismogenic zone, may be smaller than suggested by past models while also varying along-strike.