GOTHIC Survey: Dual AGN in Merging Galaxies

• GOTHIC Survey is an observational program that detects dual and multiple active galactic nuclei in merging galaxies using automated morphology-based methods.
• It employs a graph-boosted iterated hill-climbing algorithm on multi-band SDSS imaging, combined with spectroscopic filters, to isolate candidate nuclei.
• Results indicate that dual AGN primarily inhabit massive, classical bulge hosts, supporting merger-driven models of supermassive black hole growth and galaxy evolution.

Searching arXiv for recent GOTHIC survey papers and related uses of the acronym. arXiv Search Query: GOTHIC survey dual AGN Bhattacharya arXiv OR "Automated Detection of Double Nuclei Galaxies using GOTHIC" The GOTHIC survey is an observational program designed to systematically identify and characterize close multiple-nucleus galaxies in wide-field imaging, with the specific aim of isolating dual and multiple active galactic nuclei in merger environments. In this usage, GOTHIC denotes “Graph-BOosted iterated HIll-Climbing,” an image-analysis pipeline applied to Sloan Digital Sky Survey data to detect double-nucleus galaxies, followed by spectroscopic and structural analyses of the detected systems. The survey has developed from automated discovery in SDSS DR16 to host-bulge decomposition in SDSS DR18 and pPXF-based spectroscopy of dual nuclei, yielding constraints on merger demographics, bulge morphology, stellar populations, metallicity, and supermassive black hole growth (Bhattacharya et al., 2020, Nehal et al., 23 Sep 2025, Biswas et al., 18 Jun 2026).

1. Definition, scope, and nomenclature

Within extragalactic astronomy, the GOTHIC survey refers to a program for finding close pairs and higher-order multiples of galactic nuclei and for determining which of those systems host dual AGN. Its scientific motivation is to use galaxy mergers as laboratories for studying SMBH growth, AGN feedback, and bulge/galaxy co-evolution. A central question is whether dual AGN preferentially inhabit classical bulges or pseudobulges, and how the structural and spectroscopic properties of the host nuclei encode merger stage and fueling history (Nehal et al., 23 Sep 2025).

The survey’s operational distinction between double-nucleus galaxies and dual AGN is essential. A double-nucleus galaxy is a morphology-selected system with two or more closely separated nuclei in imaging. A dual AGN is a stricter spectroscopic subset in which two nuclei are classified as AGN from emission-line diagnostics. This distinction explains why the parent GOTHIC catalog is much larger than the dual-AGN subset (Bhattacharya et al., 2020).

The acronym is not unique across arXiv. “GOTHIC” also appears in unrelated contexts, including a gravitational oct-tree GPU code and work on the “gothic locus” in Teichmüller dynamics. Those usages are separate from the SDSS-based merger and dual-AGN survey discussed here (Miki et al., 2016, Torres-Teigell, 2019, Möller et al., 2018).

2. Detection algorithm and survey construction

At the core of the survey is an automated image-domain detector for multiple nuclei. The algorithm combines “graph boosting,” which maps a smoothed galaxy image to an intensity-weighted graph and accentuates local maxima, with “iterated hill-climbing,” which follows intensity gradients to convergence points identified as candidate nuclei (Bhattacharya et al., 2020). In the survey implementation, the input is a $40''\times40''$ SDSS cutout in each band. After Gaussian smoothing, each pixel becomes a graph node with initial weight

$W_0(v_p)=I(p),$

and neighbouring pixels are connected with edge weights such as

$w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$

Peak enhancement is performed through an iterative update,

$$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$

with a small learning rate $\alpha\sim0.1$. Hill-climbing is then launched from many seeds in the high-intensity region determined from a Sérsic-profile cutoff,

$$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$

and each walker moves according to

$$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$

until no higher-weight neighbour exists (Bhattacharya et al., 2020).

Candidate pruning is rule-based. The filter requires exactly two peaks within the central envelope, consistent detections across the $g,r,i,u$ bands within $\pm1$ pixel in $(x,y)$, and a separation range of approximately $W_0(v_p)=I(p),$0 to $W_0(v_p)=I(p),$1 for double-nucleus candidacy. In the original survey design, preprocessing also included Canny edge detection to isolate the galaxy envelope and a pixel-histogram Sérsic fit to estimate $W_0(v_p)=I(p),$2, $W_0(v_p)=I(p),$3, and $W_0(v_p)=I(p),$4 (Bhattacharya et al., 2020).

The initial blind run was carried out on one million SDSS DR16 galaxies with spectroscopy. In the parent search, galaxies were drawn uniformly at random from SDSS DR16 over approximately $W_0(v_p)=I(p),$5, using $W_0(v_p)=I(p),$6 imaging, with the $W_0(v_p)=I(p),$7 band dropped because of higher noise. The algorithm flagged 104,412 double-nucleus candidates, which were reduced to 949 after six band-consistency and spectroscopic filters, and manual vetting yielded 681 confirmed double-nucleus systems (Bhattacharya et al., 2020).

3. Catalog layers, selection logic, and validation

The survey pipeline is stratified: morphology first, spectroscopy second, and detailed follow-up last. This layered design is visible in the differing sample sizes reported by the discovery, bulge-structure, and spectroscopic studies. Those samples are not contradictory; they correspond to progressively stricter cuts and distinct downstream analyses (Bhattacharya et al., 2020, Nehal et al., 23 Sep 2025, Biswas et al., 18 Jun 2026).

Stage	Selection basis	Reported sample
Parent morphology catalog	Automated detection + manual vetting	681 confirmed multi-nuclei systems
Dual-AGN candidate set	BPT-selected AGN-ionized nuclei from close pairs	159 dual AGN candidates
Bulge-structure subset	Size cut, GALFIT convergence, cleaning	64 reliable systems, totaling 131 bulges
Spectroscopic quality subset	SDSS spectra with S/N and pPXF quality cuts	915 nuclei, including 341 dual systems with both components good

Algorithm validation in the discovery paper used 47 bona fide double-nucleus galaxies and 47 random single-nucleus controls. The reported confusion matrix gave 47/47 true positives, 0/47 false negatives, 7/47 false positives, and 40/47 true negatives, corresponding to 100% completeness and 85% purity on that benchmark (Bhattacharya et al., 2020).

From the 681 confirmed multi-nuclei systems, the later spectroscopic study described the composition of the full morphology-selected catalog as 652 dual, 27 triple, and 2 quadruple systems, with a redshift distribution peaking at $W_0(v_p)=I(p),$8 and SDSS $W_0(v_p)=I(p),$9 mag. After requiring continuum $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$0 and acceptable pPXF quality flags, 915 nuclei remained, including 341 dual systems in which both components satisfied the spectroscopic quality criteria (Biswas et al., 18 Jun 2026).

A common misconception is to identify all double-nucleus systems with dual AGN. The survey results do not support that equivalence. In the discovery study, 1,098 nuclei in the 681 double-nucleus galaxies had reliable line measurements, and 581 lay in AGN regions of the diagnostic plane, but only systems with two AGN-classified nuclei were retained as dual AGN, yielding 159 such systems (Bhattacharya et al., 2020). The spectroscopic work sharpened this distinction further by reporting pair classes among the 99 dual systems where both nuclei could be classified: 31 AGN–AGN, 7 SF–SF, 13 Comp–Comp, 8 AGN–SF, 34 AGN–Comp, and 6 Comp–SF (Biswas et al., 18 Jun 2026).

4. Host-bulge morphology and structural decomposition

The bulge-morphology study used SDSS DR18 calibrated FITS cutouts in the $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$1, $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$2, and $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$3 bands, all at pixel scale $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$4. The reduction protocol extracted $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$5 pixel cutouts centered on each nucleus, converted nanomaggies pixel$w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$6 to counts pixel$w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$7 using the NMGY header keyword, added a constant sky pedestal of 1000 counts for GALFIT stability, constructed a field PSF from nearby unsaturated stars, masked foreground stars and background galaxies when needed, and excluded galaxies with de Vaucouleurs radius smaller than $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$8 to ensure resolvable structure (Nehal et al., 23 Sep 2025).

Two-dimensional decomposition was performed with GALFIT. Each image was modeled with one or more Sérsic components for the bulges or nuclei, with free Sérsic index $w(u,v)=\tfrac12\bigl(I(u)+I(v)\bigr).$9, effective radius $$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$0, axis ratio, and position angle. If a disc was evident, an exponential component with Sérsic $$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$1 was added. In systems with multiple nuclei, one Sérsic bulge component was fit per nucleus. The adopted surface-brightness law was

$$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$2

with $$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$3 satisfying

$$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$4

Bulge luminosities came from GALFIT magnitudes in $$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$5, $$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$6, and $$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$7; $$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$8 colors were measured after Galactic extinction correction and K-correction; stellar mass-to-light ratios used

$$W_t(v)\;=\;W_{t-1}(v)\;+\;\alpha\sum_{u\in N(v)}w(u,v)\,\bigl[W_{t-1}(u)-W_{t-1}(v)\bigr],$$9

and bulge stellar masses were computed from

$\alpha\sim0.1$0

For ellipticals, the single-component $\alpha\sim0.1$1 was interpreted as the galaxy effective radius (Nehal et al., 23 Sep 2025).

Classification was driven primarily by Sérsic index and disc presence. Pseudobulges were defined by $\alpha\sim0.1$2, classical bulges by $\alpha\sim0.1$3, and ellipticals as single-component fits with no exponential disc and typically $\alpha\sim0.1$4 (Nehal et al., 23 Sep 2025).

The final structural sample comprised 64 reliable systems, specifically 61 dual-AGN host galaxy pairs plus 3 confirmed triplets, totaling 131 bulges at redshifts $\alpha\sim0.1$5–$\alpha\sim0.1$6. In the $\alpha\sim0.1$7 band, the Sérsic indices spanned $\alpha\sim0.1$8 with median $\alpha\sim0.1$9. Across $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$0, approximately 70–80% of bulges had $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$1, while approximately 20% had $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$2. Table 4 in that study summarized the morphological breakdown as 105 classical bulges and 26 pseudobulges; 73 occurred in single-component elliptical fits and 58 in disc galaxies (Nehal et al., 23 Sep 2025).

Bulge masses ranged from $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$3 to $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$4, with median $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$5. Elliptical hosts reached effective radii up to approximately 21 kpc, whereas bulge components in disc galaxies had $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$6 kpc. The merger morphological mix, reported for 60 dual-AGN systems excluding triplets, was 24 elliptical–elliptical systems, 24 elliptical–disc systems, and 18 disc–disc systems; the same study summarized this as roughly 70% of dual AGN involving at least one elliptical bulge host (Nehal et al., 23 Sep 2025). The numerical tension between the listed counts and the stated percentage is present in the reported data; a plausible implication is that the authors intended the qualitative conclusion—dominance of systems containing at least one elliptical—to be read more strongly than the raw tabulated arithmetic.

5. Spectroscopic characterization of dual nuclei

The spectroscopic follow-up used SDSS spectra and the penalized pixel-fitting code pPXF in a two-stage workflow. Stage A performed a stellar-only fit with gas lines masked to determine stellar kinematics. Stage B performed a full stars-plus-gas fit to derive stellar populations and emission-line fluxes. Common steps included median flux normalization, rescaling of the noise spectrum to enforce $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$7 through a wild-bootstrap procedure, the use of the E-MILES SSP library over approximately $$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$8–$$I(r)=I_e\exp\!\bigl[-b_n\bigl((r/r_e)^{1/n}-1\bigr)\bigr],\quad I(r)\ge I_{\rm cut}\equiv I_e\,e^{-b_n},$$9 Å, additive and multiplicative Legendre polynomials for continuum corrections, and robust $$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$0-clipping of residuals (Biswas et al., 18 Jun 2026).

Stellar kinematics were fit over the full SDSS range, about $$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$1–$$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$2 Å, excluding masked emission regions. The primary outputs were $$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$3 and $$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$4, with typical uncertainties $$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$5 from 50 bootstrap realizations. Emission-line measurements in the second fit included H$$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$6 4861, [O III] $$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$7, H$$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$8 6563, and [N II] $$v\;\rightarrow\;\arg\max_{u\in N(v)}W_{\rm boosted}(u),$$9, with classification requiring all four key lines at at least $g,r,i,u$0 significance (Biswas et al., 18 Jun 2026).

BPT classification followed the standard $g,r,i,u$1 versus $g,r,i,u$2 plane with the Kauffmann et al. 2003 lower boundary, the Kewley et al. 2006 upper boundary, and the Cid Fernandes et al. 2010 LINER–Seyfert divider, producing four classes: star-forming, composite, AGN–LINER, and AGN–Seyfert. Median light-weighted metallicities were reported as $g,r,i,u$3 for star-forming nuclei, $g,r,i,u$4 for composite, $g,r,i,u$5 for AGN–LINER, and $g,r,i,u$6 for AGN–Seyfert. Median light-weighted ages $g,r,i,u$7 were 8.34, 9.20, 9.86, and 9.83 for the same sequence of classes (Biswas et al., 18 Jun 2026).

Stellar masses were estimated by two methods, of which the adopted one used the $g,r,i,u$8-band apparent magnitude plus the pPXF-derived $g,r,i,u$9, with a typical $\pm1$0 uncertainty of about 0.1 dex from bootstrap realizations. SMBH masses were then derived from the $\pm1$1–$\pm1$2 relation

$\pm1$3

equivalently $\pm1$4. The Eddington luminosity was written as

$\pm1$5

Median $\pm1$6 values were 11.27, 12.10, 12.73, and 12.81 for star-forming, composite, AGN–LINER, and AGN–Seyfert nuclei, respectively (Biswas et al., 18 Jun 2026).

Several comparative trends were emphasized. For dual systems, the best-fit black-hole-mass versus stellar-mass relation had slope $\pm1$7 and intercept $\pm1$8; for AGN–Seyferts alone, slope $\pm1$9 and intercept $(x,y)$0; for AGN–LINERs, slope $(x,y)$1 and intercept $(x,y)$2; and for the $(x,y)$3 subset, slope $(x,y)$4 and intercept $(x,y)$5. The recalculated single-nucleus comparison sample from Reines and Volonteri 2015 gave slope $(x,y)$6 and intercept $(x,y)$7. The reported result was that, at fixed $(x,y)$8, dual AGN lie systematically above the single-nucleus relation, implying black-hole masses approximately twice as high. The paper also noted 169 dual systems with $(x,y)$9 despite $W_0(v_p)=I(p),$00, and found no clear dependence of $W_0(v_p)=I(p),$01 on projected separation in the 5–35 kpc range or on major/minor merger status (Biswas et al., 18 Jun 2026).

6. Interpretation, scientific significance, and limitations

Across the survey sequence, a consistent empirical picture emerges. In the discovery paper, color–magnitude trends placed SF–SF pairs in the blue cloud with median $W_0(v_p)=I(p),$02, AGN–SF pairs in the green valley, and dual AGN on the red sequence with $W_0(v_p)=I(p),$03; as nucleus separation shrank, hosts reddened by $W_0(v_p)=I(p),$04 mag. The same study reported that AGN are ubiquitous in systems with $W_0(v_p)=I(p),$05 kpc and that the demographic rate of dual AGN was approximately $W_0(v_p)=I(p),$06 of all galaxies, far below naive expectations from merger rates, suggesting either a short duty cycle for simultaneous AGN activity or substantial obscuration and host dilution (Bhattacharya et al., 2020).

The structural study sharpened the host-galaxy side of that picture by showing that dual AGN overwhelmingly inhabit massive classical bulges or ellipticals, rather than pseudobulges, and by associating the systems primarily with evolved, red, quenched hosts and major mergers. It reported that major mergers with galaxy mass ratio $W_0(v_p)=I(p),$07 account for approximately 60% of the systems and interpreted the prevalence of classical bulges and ellipticals as evidence that dual SMBH accretion episodes are preferentially triggered in merger-built stellar potentials rather than in secularly evolved pseudobulges (Nehal et al., 23 Sep 2025).

The spectroscopic work extended the interpretation to black-hole growth and chemical evolution. Its stated conclusion was that SMBH masses are higher for black holes in galaxy mergers than for single nuclei at fixed stellar mass, revealing that SMBHs grow during the galaxy merging process and not only due to the merger of SMBHs. It also reported mass–metallicity offsets relative to Gallazzi et al. 2005: low-mass dual nuclei with $W_0(v_p)=I(p),$08 had median $W_0(v_p)=I(p),$09 versus $W_0(v_p)=I(p),$10 in Gallazzi et al., while high-mass systems with $W_0(v_p)=I(p),$11 had median $W_0(v_p)=I(p),$12 versus $W_0(v_p)=I(p),$13. The interpretation offered was that mergers induce mixing and inflows that alter central metallicities (Biswas et al., 18 Jun 2026).

Several limitations are explicit. The bulge study emphasizes that its sample is limited, with only 69 systems converging in GALFIT and 64 retained after further cleaning (Nehal et al., 23 Sep 2025). The spectroscopic analysis is restricted to SDSS spectra satisfying continuum $W_0(v_p)=I(p),$14, and robust BPT classification requires all four key lines at at least $W_0(v_p)=I(p),$15 significance (Biswas et al., 18 Jun 2026). More broadly, the pipeline’s morphology-first logic means that the survey is optimized for close, resolvable multiple nuclei in SDSS imaging rather than for every physically dual SMBH system.

Taken together, the GOTHIC survey defines a coherent observational program: automated morphology-based discovery of close multiple nuclei, spectroscopic isolation of dual AGN, and quantitative characterization of host structure and nuclear properties. The combined results support a merger-driven scenario in which simultaneous SMBH fueling, bulge assembly, chemical redistribution, and star-formation quenching are linked phases of hierarchical galaxy evolution (Bhattacharya et al., 2020, Nehal et al., 23 Sep 2025, Biswas et al., 18 Jun 2026).