EconBase
← Back to paper

Making Interpretable Discoveries from Unstructured Data: A High-Dimensional Multiple Hypothesis Testing Approach

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

117,076 characters · 18 sections · 147 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Making Interpretable Discoveries from Unstructured Data: A High-Dimensional Multiple Hypothesis Testing Approach

abstractSocial scientists are increasingly turning to unstructured datasets to unlock new empirical insights, e.g., estimating descriptive statistics of or causal effects on quantitative measures derived from text, audio, or video data. In many such settings, unsupervised analysis is of primary interest, in that the researcher does not want to (or cannot) manually pre-specify all important aspects of the unstructured data to measure; they are interested in “discovery.” This paper proposes a general and flexible framework for pursuing such discovery from unstructured data in a statistically principled way. The framework leverages recent methods from the literature on AI interpretability to map unstructured data points to high-dimensional, sparse, and interpretable “concept embeddings”; computes statistics from these concept embeddings for testing interpretable, concept-by-concept hypotheses; performs selective inference on these hypotheses using algorithms validated by new results in high-dimensional central limit theory, producing a selected set (“discoveries”); and both generates and evaluates human-interpretable natural language descriptions of these discoveries. The proposed framework has few researcher degrees of freedom, is robust to data snooping and other post-selection inference concerns, and facilitates fast and inexpensive sensitivity analysis and replication. Applications to recent descriptive and causal analyses of unstructured data in empirical economics are explored. Open source code is provided for researchers to implement the framework in their own projects.

Introduction

Empowered by recent developments in machine learning and AI, researchers in the social sciences are increasingly leveraging sources of unstructured data---such as text, audio, images, and videos---in quantitative analyses. In economics, interest in unstructured data sources is especially widespread, ranging from using speech recordings from FOMC meetings to better understand monetary policy gorodnichenko_voice_2023; to using videos of start-up pitches to study entrepreneurship and investment hu_persuading_2025; to using open-ended survey questions to probe the economic behavior and beliefs of individuals haaland_understanding_2024; to using mug shots to study judicial decisionmaking ludwig_machine_2024; to using visual art or written narratives to infer long-run living standards gorin_state_2025,lagakos_american_2025; to using qualitative interview transcripts to better understand the impacts of RCT treatments bergman_creating_2024,krause_impact_2025.

Unstructured data afford new opportunities for measuring social phenomena that would be otherwise unmeasurable. In some settings, theory or prior empirical research suggests important concepts that unstructured data newly facilitates measures of; in others, the promise of unstructured data is that it allows for the discovery of entirely novel concepts (and associated measures) with social or economic relevance (c.f., ludwig_machine_2024). When modern AI/ML methods are brought to bear on these goals, the former typically falls under the heading of “supervised learning,” whereby the researcher makes predictions of their ex-ante known quantity of interest using the unstructured data as features, and the latter typically falls under the heading of “unsupervised learning,” which is less amenable to a simple prediction-based framework, and demands that AI/ML methods uncover latent structure in the unstructured data with minimal researcher input. Though new statistical and econometric frameworks have recently become available for interpretable and statistically principled analyses of unstructured data in the supervised learning setting (e.g., angelopoulos_prediction-powered_2023,ludwig_large_2024,carlson_unifying_2025), it is an open question as to how one best performs interpretable and statistically principled analyses of unstructured data in unsupervised learning settings emphasizing discovery. The framework proposed in this paper provides one possible answer to this open question.

Consider two motivating examples for the proposed framework, one focused on descriptive analysis, and another focused on causal inference. Both consider text as the unstructured data modality of interest, though the proposed framework applies more generally to other data modalities (see, e.g., the AI interpretability literature review in Section (ref)).

example[Information treatment RCT, e.g., bursztyn_justifying_2023] A researcher is interested in running a randomized controlled trial (RCT) with an information treatment. A key outcome of the experiment is the participant's subjective reaction to the information treatment, and the researcher would like to discover any systematic differences in reactions across the treatment and control groups. As such, for all participants, the researcher collects text responses to open-ended questions eliciting reactions to the experimental intervention, and is interested in learning any relevant causal effects on the presence of concepts surfaced in these text responses.
example[Open-ended survey analysis, e.g., stantcheva_why_2024] A researcher is interested in understanding the attitudes of some population of interest towards some particular (economic) policy. To learn about population attitudes without priming respondents to discuss particular topics, the researcher runs a survey on a representative sample from the population of interest that asks many open-ended questions about attitudes and beliefs, to which survey respondents reply in text. The researcher is interested in describing any frequently discussed concepts that are surfaced from this open-ended survey.

To address these examples---among many others---this paper proposes a general framework for conducting interpretable and rigorous statistical inference on unstructured data per the following four steps:

enumerate• Creating concept embeddings: The researcher begins with an unstructured dataset they are interested in analyzing. Each unstructured data point is mapped to a “concept embedding,” which is a binary vector.\footnote{The name “concept embedding” is chosen to be consistent with other parts of the CS and AI interpretability literatures that make use of similar objects, e.g., bhalla_interpreting_2024,jiang_interpretable_2025.} Each element of the vector indicates the presence (or lack thereof) of a single human-interpretable feature---or, informally, a “concept”---in that unstructured data point. The creation of these interpretable concept embeddings can be accomplished using pretrained dictionary learning models newly popular in the AI interpretability literature for deep neural networks (DNNs) such as large language models (LLMs). These concept embedding vectors are high-dimensional, and catalog a vast number of concepts. Heuristically, a concept embedding implemented with a dictionary learning model should be thought of as an inventory of all concepts that, e.g., a LLM “needed to learn” for the purpose of good next token prediction over a massive internet-scale corpus. By construction, these concept embedding vectors are also sparse, i.e., only the few most “meaningful” concepts represented in the unstructured data point are cataloged as present by the dictionary.\footnote{Meaningfulness, in a technical sense, relates to the reconstruction loss associated with learning the dictionary feature under a sparsity penalty; for more on this, see, e.g., bricken2023monosemanticity.} \begin{examplecont}{ex:rct} Each experimental participant's post-treatment text response is passed through a LLM equipped with pretrained sparse autoencoder (SAE), a popular implementation of a dictionary learning model bricken2023monosemanticity. The SAE's feature activations are pooled and thresholded at zero to serve as a concept embedding. There is one concept embedding of each text response (and one text response per experimental participant), each of which is of dimension on the order $p\approx 10^5$, though there are only on the order of $n\approx 10^2$ participants across both treatment and control groups. Each concept in the embedding can be given an English language description through a process called “automatic interpretation” (discussed in step 4), e.g., concept number 1 in the concept embedding may indicate “the presence of sentences ending with exclamation points,” concept number 2 may indicate “the discussion of politics,” and so on.\footnote{For more examples of feature descriptions for, e.g., the Gemma Scope family of SAEs lieberum_gemma_2024, see: https://www.neuronpedia.org/gemma-scope.} \end{examplecont} • Formulating concept-level hypotheses: The researcher, now with access to a new dataset of concept embeddings, formulates $p$ concept-specific hypothesis tests and associated test statistics---one for each concept cataloged in the concept embedding. The researcher has great flexibility in designing relevant, concept-level hypotheses and t-stats at this stage of the framework, which may employ data sources beyond the concept embeddings themselves. \begin{examplecont}{ex:rct} The researcher computes a difference in means across treatment and control groups for each concept indicator (a binary entry in the concept embedding) to test a total of $p$ null hypotheses of zero average treatment effect (ATE) on the probability of a given concept appearing in a text response. \end{examplecont} • \textit{High-dimensional multiple hypothesis testing}: All $p$ concept-specific hypotheses are tested using statistical procedures with high-dimensional selective error guarantees validated by novel results in high-dimensional multiple hypothesis testing and central limit theory. These results prove that the algorithms of romano_control_2007 for controlling the $k$ familywise error rate ($k$-FWER)---the probability of making $k$ or more false rejections---are valid in high-dimensions under mild conditions for small $k$. $k$-FWER control for $k>1$ is less conservative than simple FWER control---control over the probability of making \textit{any} false rejections\footnote{Note that FWER is equivalent to $k$-FWER when $k=1$.}---leading to more powerful tests that permit the researcher to make more discoveries (rejections of true nulls). Moreover, the researcher may invert these tests to form generalized simultaneous confidence intervals, i.e., the probability that $k$ or more parameters fall outside the simultaneous intervals is less than or equal to some desired level. \begin{examplecont}{ex:rct} An algorithm based on the work of romano_control_2007 is deployed to test each of the $p$ null hypotheses of no average treatment effect, guaranteeing $k$-FWER control at level $\alpha$ (in large samples) for a researcher chosen value of $k$ (that is sufficiently small), even though $p \gg n$. Bootstrapped critical values may also be used to form generalized simultaneous confidence intervals for the difference-in-means point estimates of the average treatment effects. \end{examplecont} • \textit{Generating and evaluating automatic interpretations}: The output of the framework so far is a selected set of rejected, concept-level null hypotheses (and corresponding generalized simultaneous confidence intervals and point estimates), a.k.a., “discoveries.” By virtue of the dictionary learning technology used to create them, however, the concepts encoded by the concept embedding are only indexed numerically, and do not come with natural language descriptions at first---we know, e.g., the first entry of the concept embedding likely encodes a human-interpretable concept, and that we rejected the null hypothesis associated with it, but we do not yet have an English description for what that concept is. As such, in the final stage of the framework, the researcher generates refined natural language descriptions of the concepts in the concept embedding using a variant of a popular \textit{automatic interpretation} (or “autointerp”) method that is tailored to the unstructured dataset being explored.\footnote{Many pretrained dictionary learning models come with “off-the-shelf” feature descriptions, but they are usually of poor quality, and benefit from refinement based on information in the dataset of interest; this is discussed further in Section (ref).} In order to rigorously evaluate the quality of these autointerp descriptions, we propose a new statistical formalization of a workhorse autointerp evaluation method from paulo_automatically_2024, along with new estimators that permit statistically and conceptually principled autointerp quality evaluation. \begin{examplecont}{ex:rct} The researcher now has in hand some subset of the $p$ null hypotheses that have been rejected: a set of indices associated with particular elements of the concept embeddings. In order to provide a natural language description for each of the selected concepts for which there is evidence of a nonzero average treatment effect (which are otherwise just identified by their numeric index), a dedicated “explainer LLM” (which may be closed or open source) reads examples of texts in which a given feature is present (or presents with high SAE activations) and generates an English interpretation for what concept is encoded. The credibility of this LLM-derived description of the concept is then evaluated with a specially designed classification task on a held-out evaluation sample using newly proposed estimators. \end{examplecont}

The final output of the proposed framework is therefore: a set of concept-level discoveries (rejections of nulls) that come with generalized familywise error rate guarantees; natural (e.g., English) language descriptions of those discoveries; and interpretable scores measuring the quality of those descriptions. A schematic for the framework is depicted in Figure (ref).

examplecont{ex:rct} The researcher discovers that the null of no average treatment effect is rejected for concepts indexed $\{4, 101, 5030, ...\}$. Concept 4 is described by an explainer LLM as activating on “mentions of politics” and this description is estimated as having a high “quality” score, i.e., the information treatment credibly causes increased discussion of politics; and so on.
figure[figure omitted — 194 chars of source]

The status quo for many researchers interested in high quality, interpretable discovery from unstructured data typically involves an iterative and discretionary “coding” process, compromised of manually perusing the unstructured data to glean concepts or topics that appear salient (perhaps in collaboration with LLMs) and then codifying them in a “codebook” to operationalize quantitative measurement.\footnote{For just a few very recent examples of this approach in “top five” economics journals, see, e.g., andre_subjective_2022,geiecke_conversations_2024.} (See, e.g., haaland_understanding_2024 for a more detailed characterization of this process in empirical economics research.) Though intuitive, such an approach to discovery suffers from several distinct scientific and statistical problems, all of which the proposed framework aims to address:

enumerate• The streetlight effect: The “streetlight effect” describes the cognitive bias that occurs when researchers only think to measure that which is easiest or most obvious to measure (idiomatically speaking, the researcher is “only searching for their lost keys where the light is shining.”\footnote{The same problem is sometimes cast as failing to measure or discover important “unknown unknowns.”}) Even if a researcher is well-intentioned, they may simply fail to notice scientifically important patterns that appear in unstructured datasets when manually inspecting unstructured data. In the proposed framework, AI interpretability methods are repurposed to measure the presence of a vast number of human-interpretable concepts in an unstructured dataset---every concept or feature that a modern LLM needed to learn to excel at pretraining on an internet-scale corpus---reducing the possibility that there are any truly important measures of concepts ignored by the researcher (and thus, in a sense, shining a much more powerful streetlight on the data). • Researcher degrees of freedom: Because of the incentives researchers face, any largely opaque, ad hoc, or discretionary process of defining concepts for later measurement runs the risk of being contaminated by motivated data mining, e.g., only defining and reporting measures that support a preordained conclusion. As such, any procedure for conducting discovery from unstructured data with too many researcher degrees of freedom suffers from less credibility than one that ties the researcher's hands in some transparent way. The proposed framework is automatic in execution, and does not require a “human-in-the-loop” to direct what is and is not measured: in other words, it has very few researcher degrees of freedom. The few choices a researcher must make in the proposed framework, such as a choice of dictionary learning model and explainer LLM, can also easily be preregistered, further hampering the possibility that a researcher applies this framework with many different models to select on favorable results. Even without preregistration, cherry-picking dictionary learning models to reverse engineer a desired result is unlikely to succeed given that, e.g., SAE features tend to exhibit some degree of universality across architectures (c.f., wang2025towards). • Data snooping and post-selection inference: Even when data snooping is not motivated by some desired conclusion, both defining and taking measurements based on the same observed dataset realization leads to well known problems in post-selection inference, which invalidates downstream statistical conclusions if handled incorrectly.\footnote{See, e.g., kuchibhotla_post-selection_2022 for more on the problems of and solutions to post-selection inference.} The proposed framework relies on using data-independent, ex-ante curation of a high-dimensional set of meaningful concept measurements---coupled with appropriate high-dimensional selective inference procedures---to achieve many of the scientific benefits of exploration without running into the post-selection inference problems introduced by data snooping. • Replicability: Even for thoughtful discretionary exploratory analyses of unstructured data, replication is often costly (or infeasible, without proper documentation), and assessing sensitivity to alternative analysis choices is burdensome. By contrast, the automaticity of the proposed framework makes it easy to reproduce discoveries and their descriptions and analyze the sensitivity of results with respect to those few choices a researcher must make (e.g., choice of concept embedding, or explainer LLM). Further facilitating replication and sensitivity analysis is that the proposed framework is inexpensive both in terms of compute cost and time investment: both empirical examples in this paper were computed on Google Colab notebooks using a single A100 GPU, access to which only costs \$10 per month (and is free for students and academic faculty), and takes on the order of minutes to run using this hardware.\footnote{As of April 2026.}

Though the use of sample splitting in the process of discovery from unstructured data (as advocated for in, e.g., modarressi_causal_2025) addresses the important post-selection inference problem, it does not directly address the other problems discussed above. Moreover, sample splitting reduces the statistical power available for rejecting nulls, which may already be precious in the small (especially experimental) datasets common to the social sciences.\footnote{Note that, in most applications concerning discovery from unstructured data that involve human discretion over choice of measurements, efficiency cannot be recovered from cross-fitting.} The proposed framework also incurs a reduction in power due to multiple hypothesis testing corrections, though such corrections would still be necessary under sample splitting if many hypotheses were transported across samples for the purposes of evaluation---as may often be of interest in settings emphasizing discovery.

It is worth noting that problems surrounding post-selection inference, researcher degrees of freedom, and motivated data mining have been taken very seriously in other literatures, such as the analysis of RCTs in economics, where these concerns motivated the creation of preanalysis plans (PAPs) for experiments. One analogy for the proposed framework, then, is that it provides similar value to a PAP for the setting of discovery from unstructured data. Just like a PAP, the proposed framework can be accompanied by contingent, complementary analyses of the unstructured dataset, so long as these analyses are evaluated by readers with the appropriate skepticism banerjee_praise_2020. However, unlike writing a comprehensive PAP for an experiment, the proposed framework places very little burden on the researcher: the same default template for analysis (the concepts measured by the SAE) can be reused across many diverse empirical settings; in fact, the credibility of analysis is bolstered by not deviating from the default choice of SAE.\footnote{When the unstructured data being analyzed literally comes from an experiment, preregistration of the proposed framework in an actual PAP is also easily declared (it simply involves stipulating a particular choice of LLM and SAE under the proposed concept embedding).}

Though useful for experiments, it is also important to note that a comparative advantage of the proposed framework may lie in its application to observational datasets. In experimental settings (especially when experiments are cheap to run), it may not be necessary to pursue statistically principled discovery by carefully leveraging two stages of analysis, following the blueprint of ludwig_machine_2024: the first stage can perform scientifically informative but not necessarily statistically valid discovery of important concepts (e.g., in a pilot experiment), and the second stage can implement a more powerful experiment that formally assesses a few important concepts surfaced from the first stage. However, in settings featuring observational data, or even certain experimental settings, the collection of a new, independent dataset to be used solely for formal evaluation is often infeasible (e.g., a historical setting for which only one sample exists, or a very expensive large-scale experiment), making statistically valid exploration within a single dataset crucial---something accomplished by the proposed framework. That said, even when viable, two stage procedures for discovery and testing suffer from some of the same scientific problems as sample splitting on a single dataset, owing to discretion over what concepts are ported into the second stage for measurement. The proposed framework therefore still adds value to readers concerned with cherry-picking and motivated data mining across stages. Furthermore, the question of how to structure open-ended discovery in the first stage---even if not statistically grounded---is an open one, and the proposed framework provides support for the use of pretrained SAEs to discover concepts (complementing work by movva_sparse_2025,wang_your_2026 on the value of learning SAEs for concept discovery).

A core intended contribution of this paper is its development of a general-purpose framework for making statistically and scientifically principled discoveries from unstructured data, offering solutions to many problems with status quo practice. However, in order to operationalize this framework, we also make the following theoretical contributions, which themselves may be of independent interest:

enumerate• Extensions of two frontier econometric toolkits: extending the high-dimensional “many approximate means” framework of belloni_high-dimensional_2018 to the study of the $k$-th largest coordinate of sums of independent random vectors, and extending the application of the $k$-FWER controlling procedures of romano_control_2007 to the high-dimensional setting. • A new conceptual and statistical formalization of a popular framework for autointerp scoring from paulo_automatically_2024, including novel estimation procedures.

The remainder of this paper is structured as follows: Section 2 discusses related literature; Section 3 describes the framework in detail; Section 4 applies the framework to two recent papers making use of unstructured data for discovery bursztyn_justifying_2023,stantcheva_why_2024; and Section 5 concludes.

Related Literature

The framework proposed in this paper is related to recent works from literatures spanning economics, statistics, and computer science.

\paragraph{Econometric methods for unstructured data.} Motivated by the capabilities of modern AI/ML methods for learning from unstructured economic data (see, e.g., mullainathan_machine_2017,gentzkow_text_2019,ash_text_2023,dell_deep_2025), new econometric and statistical frameworks have been developed to facilitate principled statistical inference on low-dimensional features (predictions) learned from unstructured datasets under supervision, e.g., angelopoulos_prediction-powered_2023,ludwig_large_2024,carlson_unifying_2025,rambachan_program_2024. This most recent wave of econometric literature emphasizes nonparametric frameworks compatible with black-box AI models and that center debiasing methods, as opposed to model-based approaches, e.g., gentzkow_measuring_2019,battaglia_inference_2024. Though the framework presented in this paper is concerned with principled inference on unstructured data without making parametric assumptions, it differs from this existing econometric literature in that its goal is primarily unsupervised discovery, as opposed to supervised detection, using AI/ML models.

The proposed framework is most similar to existing works such as ludwig_machine-learning_2017,modarressi_causal_2025, both of which provide formal statistical frameworks for discovery from high-dimensional economic data. Both of these works focus specifically on causal inference in RCTs, unlike the proposed framework, which generally treats any inference question that can be posed as a hypothesis test (including causal inference based on RCT data). Both ludwig_machine-learning_2017 and modarressi_causal_2025 employ the intuition of a “reverse regression” to formally test whether a distribution of (high-dimensional) outcomes is the same across treatment and control groups as a means of discovery. However, both papers note that such a coarse test is often insufficient for developing scientific insight into the causal effects of treatment, motivating methods for probing what changed across the outcome distributions. modarressi_causal_2025 specifically explores this idea for (unstructured) text data, suggesting a framework that splits the researcher's sample: in one split, a human-LLM hybrid pipeline is used to hypothesize a low-dimensional set of “causal themes” that differ across treatment and control, and in the other split a human researcher scores texts according to these themes and evaluates differences in means. As such, for the setting of inference in RCTs, the framework of modarressi_causal_2025 differs from the proposed framework along several dimensions, including: encouraging human intervention to guide the measurement of themes/concepts (as opposed to minimizing the influence of a human-in-the-loop), using LLM reasoning as a low-dimensional measurement tool (as opposed to leveraging AI interpretability methods like SAEs for high-dimensional measurement), and invoking sample splitting for valid inference (instead of relying on many independently curated concepts and multiple hypothesis testing corrections).\footnote{Given the amount of human intervention involved, modarressi_causal_2025 discuss concerns pertaining to the replicability of and researcher degrees of freedom in their framework in Section 8 of their paper.}

\paragraph{High-dimensional selective inference.} There is an extensive literature on multiple hypothesis testing (or “selective inference” more broadly) in statistics, biostatistics, and economics (see, e.g., romano_hypothesis_2010 for a review in econometrics). The methods in this literature span asymptotic and finite sample valid frameworks, and low- and high-dimensional settings. In particular, for the problem considered in this paper, the literature on asymptotically valid, high-dimensional selective inference is most relevant, as such methods permit making very few assumptions about the data (c.f., assuming independence or PRDS of p-values in the case of benjamini_controlling_1995) and readily incorporate resampling methods that improve power (relative to guarding against worst-case dependence, c.f., Bonferroni corrections or the method of benjamini_control_2001). The new results for high-dimensional central limit theory and multiple hypothesis testing appearing in this paper are derived by combining and extending recent work on high-dimensional central limit theory and Gaussian multiplier bootstraps for order statistics belloni_high-dimensional_2018,chernozhuokov_improved_2022,ding_gaussian_2025. The formal results on high-dimensional central limit theory in this paper are best viewed as extensions of some core results of the “many approximate means” framework of belloni_high-dimensional_2018, moving the framework from the study of the maximal order statistic to the study of the largest order statistics by leveraging new results from ding_gaussian_2025. These high-dimensional CLTs in turn support new high-dimensional ($p\gg n$) analysis of the $k$-FWER controlling step-down procedures of romano_control_2007,romano_formalized_2008---which were originally studied in a low-dimensional (fixed $p$) setting---leading to theorems that show such procedures remain valid for asymptotic $k$-FWER control in high-dimensions under several important technical conditions (chiefly, that $k$ is fixed, or “sufficiently small”). These results in asymptotic high-dimensional $k$-FWER control complement existing methods for asymptotic high-dimensional control of the FWER belloni_high-dimensional_2018 and asymptotic high-dimensional control of the FDR under weak dependence liu_phase_2014,belloni_high-dimensional_2018.

\paragraph{Hypothesis generation.} There is a growing literature in economics and computer science on hypothesis generation at a “pre-scientific” stage of empirical inquiry ludwig_machine_2024. From the perspective of this literature, the proposed framework uses dictionary learning models as a high-dimensional hypothesis generator, as similarly advocated for by movva_sparse_2025,peng_use_2025.\footnote{The framework of movva_sparse_2025 and the framework proposed in this paper differ somewhat in how they approach hypothesis generation with SAEs. movva_sparse_2025 introduces a downstream quantity of interest and discovers (selects) concepts by how well they predict this quantity under hyperparameters specified by the researcher, whereas the proposed framework uses novel methods for high-dimensional multiple hypothesis testing to discover relevant concepts. movva_sparse_2025 also emphasizes learning SAE features from dense embeddings of the researcher's dataset, whereas the proposed framework is designed to leverage pretrained, off-the-shelf SAEs that require no further training, and is thus compatible with datasets of only a few hundred data points (common in the social sciences).} However, unlike other works in this literature that suggest a second “scientific” stage of data collection and analysis to formally evaluate the fruits of statistically informal hypothesis generation, the proposed framework considers statistically principled hypothesis generation and hypothesis testing in “one shot,” given access to only one (potentially small) unstructured dataset.

\paragraph{AI interpretability.} There is a large literature in computer science on machine learning interpretability (see, e.g., doshi-velez_towards_2017). The machine learning interpretability methods leveraged by the present framework originate from a nascent though highly active literature known as “mechanistic interpretability,” which seeks to develop methods for interpreting LLM behavior via quantitative analyses of model internals (e.g., activations, weights). In particular, the dictionary learning methods for LLM interpretability implemented in this framework are called “sparse autoencoders” (SAEs), and were developed in bricken2023monosemanticity,templeton2024scaling.\footnote{SAEs are very wide autoencoders inserted at various layers of an LLM that are trained to intercept and reconstruct activations from model internals (such as the residual stream, MLPs or the like) under a sparsity-inducing penalization on a reconstruction loss. Under the “linear representation hypothesis” and “superposition hypothesis” (see, e.g., bricken2023monosemanticity for more discussion), these sparse autoencoders are thought to act as an overcomplete basis of the space of concepts in text leveraged by a LLM to make next token predictions, encouraging learned autoencoder latents to have monosemantic interpretations.} \!\footnote{There has been much debate about the success of SAEs in pursuit of AI interpretability (e.g., leask_sparse_2025), though most criticisms still support the notion that SAEs and other dictionary learning methods have a comparative advantage in discovering as opposed to detecting concepts of interest in text peng_use_2025---which is how they are applied in this framework.} Though SAEs are a promising default choice for creating concept embeddings peng_use_2025,jiang_interpretable_2025, other dictionary learning methods that have been proposed as competitors to SAEs (e.g., transcoders paulo_transcoders_2025) may similarly serve as concept embedding models, and are similarly compatible with the statistical framework proposed in this paper. The proposed framework is also related to recent independent efforts in “data-centric interpretability,” such as jiang_interpretable_2025, which leverages SAE activations to create interpretable concept embeddings of text in order to facilitate descriptive analysis of datasets relevant to LLM behavior.

SAEs and similar methods learn quantitative measures of interpretable features, but not on their own generate natural language (e.g., English) descriptions for these learned features. As such, a related literature on “automatic interpretability” (or “autointerp”) methods has become active, which seeks to coherently use LLMs to describe the features discovered by dictionary learning methods at scale. Recent important papers in this literature include: bills2023language,shaham_multimodal_2024,paulo_automatically_2024,rajamanoharan_jumping_2024. Specifically, the autointerp framework proposed in this paper adapts and statistically formalizes insights from the “detection scoring” method of paulo_automatically_2024.

Though the empirical examples considered in this paper handle one of the most popular unstructured data types---text---for which dictionary learning-based interpretability methods are the most mature, dictionary learning techniques have been successfully applied to many other modalities, including audio and images abdulaal_x-ray_2024,bhalla_interpreting_2024,fry2024towards,daujotas2024case,pluth_sparse_2025. Future work would seek to apply the proposed framework to unstructured datasets in these other data modalities.

Framework

We now consider each step of the proposed framework in detail. In what follows, we define $[p] := \{1, \dots, p\}$, and let $x \mapsto \log x$ be the natural logarithm. Further, for any vector $s \in \mathbb{R}^p$ and a subset of indices $K \subseteq[p]$, let $s_K \in \mathbb{R}^{|K|}$ denote the subvector indexed by $K$, such that $\left(s_K\right)_q=s_{j_q}$ for $j_q \in K$ and $j_1<j_2< \cdots<j_{|K|}$. We use $s \mapsto s_{[k]}$ as notation for the function that selects the $k$-th largest coordinate of the vector $s$.

Creating Concept Embeddings

Consider that the researcher has access to a dataset of size $n$, $\{(W_i, Z_i)\}_{i=1}^n$, which is sampled i.i.d. from some population of interest represented by distribution $P$.\footnote{That is, the fundamental source of uncertainty in this framework is sampling uncertainty, driven by an underlying sampling experiment in which units are sampled i.i.d. from $P$, and each unit is associated with a measurement $(W_i, Z_i)$.} The $Z_i \in \mathcal{Z}$ are unstructured data points and the $W_i \in \mathcal{W}$ are any other observed covariates of interest. The space $\mathcal{Z}$ is typically high-dimensional and semantically lean, e.g., if each $Z_i$ was a $244 \times 244$ pixel image, $\mathcal{Z}$ might then be the space of all $3 \times 244 \times 244$ arrays of RGB values, or if $Z_i$ was text that was truncated up to some maximum length, $\mathcal{Z}$ might be the space of all binary matrices of a certain dimensionality, which are concatenated one-hot encodings for each word, term, or token in the text with respect to a specific pre-defined vocabulary.

examplecont{ex:rct} The $Z_i$ is a text response elicited for experimental participant $i$; the $W_i \in \mathcal{W} = \{0,1\}$ is a treatment status indicator for experimental participant $i$.

The first stage of the framework relies on the researcher having access to a function which creates semantically rich and human-interpretable concept embeddings from unstructured data points, which we may denote $\texttt{Dict}: \mathcal{Z} \to \{0,1\}^p$. We suggest such a function can be computed by passing an unstructured data point to a deep neural network (DNN) equipped with a dictionary learning model (e.g., a LLM equipped with pretrained SAE when $Z_i$ are texts); recording the feature activations from the dictionary model; and then pooling and binarizing these feature activations to form a vector of indicators for each feature in the dictionary, indexed $j \in [p]$. There are many possible ways for a researcher to implement a $\texttt{Dict}$ function, though an especially straightforward choice for the use case of text data is implementing a $\texttt{Dict}$ function based on the SAE activations at a single layer of the DNN, for which the $j$-th entry of the output is equal to 1 if the $j$-th SAE feature activated on any token of the input text, and 0 otherwise: \[ \mathtt{Dict}(z; l)_j :=

cases&1 if SAE feature $j$ at layer $l$ activates on any token of z, \\ &0 otherwise.

\] Concrete implementation choices for the SAE architecture and layer $l$ will be discussed in the empirical examples in Section (ref).

We denote $Y_i := \mathtt{Dict}(Z_i) \in \{0,1\}^p$ (suppressing information about layer $l$ from the notation going forward) as the “concept embedding” of an unstructured data point. Each “concept indicator” $Y_{ij} \in \{0,1\}$ has the natural interpretation that a particular concept $j$ appeared in the unstructured data point $i$.\footnote{SAE feature activations are typically positively valued scalars, where magnitude is thought to correspond to some notion of intensity of activation. Other works in mechanistic interpretability consider using these raw activation values directly, perhaps with max or average pooling across tokens. However, presently, magnitude of activation is not a well-understood or highly interpretable quantitative property of dictionary learning methods, and as such the default $\texttt{Dict}$ function advocated for in the proposed framework does not incorporate this information. That said, leveraging feature activation intensity is empirically shown to be useful for deriving good autointerpretation generation strategies; see Section (ref) for more on this.}

examplecont{ex:rct} Each text response $Z_i$ is passed through an open-source LLM equipped with a pretrained SAE. The SAE activations for each token of $Z_i$ at the $l$-th layer of the LLM are max-pooled (i.e., the maximum token activation is used as the activation summary for the entire text) and then thresholded at zero, yielding a $p\approx 10^5$-dimensional binary vector $Y_i$ for each experimental participant $i$. Each $Y_{ij} \in \{0,1\}$ is a binary indicator of a particular concept being expressed in experimental participant $i$'s text response, e.g., $Y_{i1}$ could indicate whether or not experimental participant $i$'s text response contains phrases indicating uncertainty.

The $\texttt{Dict}$ function can be thought of as effectively producing a new, semantically rich and interpretable dataset of concept embeddings and associated covariates $\{(Y_i,W_i)\}_{i=1}^n$ derived from the original dataset $\{(Z_i,W_i)\}_{i=1}^n$. The crux of the proposed framework is that inference on parameters derived from each concept indicator $Y_{ij}$ (or each $(Y_{ij}, W_i)$ pair) is desirable because:

enumerate• Each $Y_{ij}$ is engineered to be “monosemantic” (indicates the presence of a single, human-interpretable concept bricken2023monosemanticity,templeton2024scaling,movva_sparse_2025,wang_your_2026) and therefore parameters derived from $Y_{ij}$ for any given $j$ (and associated hypotheses) are interpretable, e.g., $E_P[Y_{ij}]$ is the probability of human-interpretable concept $j$ appearing in an unstructured data observation in population $P$.\footnote{For more discussion of what “monosemanticity” means and its importance for LLM interpretability, see, e.g., the discussion in bricken2023monosemanticity.} • Collectively, the $p$ indicators in $Y_{i}$ for any given $i$ catalog a vast space of concepts: again, heuristically, they encode the set of all concepts that a DNN (e.g., a LLM) needed to “understand” to excel at pretraining (e.g., next token prediction) on a massive-scale corpus of unstructured data (e.g., internet-scale text corpora). As such, simultaneous inference on parameters derived from all of these concept indicators greatly reduces the possibility that a researcher misses out on discovering important “unknown unknowns” in their unstructured data; the fact that these concept measures were curated in a data-independent way further shields inference from the important scientific and statistical problems associated with data snooping.

Towards exploiting this wealth of interpretable information, we consider formulating formal hypothesis tests of parameters (estimands) of interest in the next step of the framework.

Formulating Concept-Level Hypotheses

With an interpretable dataset of concept embeddings and associated covariates $\{(Y_i,W_i)\}_{i=1}^n$, we may consider defining $p$ “one-sided” concept-level hypotheses $\{\tilde H_{0,j}\}_{j \in [p]}$, where \[\tilde H_{0,j}: \theta_j(P) := E_P[X_{ij}] \leq 0, \quad X_{ij} = h(W_i, Y_{ij}),\] for some measurable function $h$, or, alternatively, a set of “two-sided” hypotheses of interest $\{ H_{0,j}\}_{j \in [p]}$ as \[ H_{0,j}: \theta_j(P) := E_P[X_{ij}] = 0.\] As shorthand, we will denote $X := \{X_{i}\}_{i=1}^n$. Letting the concept-level parameter of interest $\theta_j(P)$ be an expectation of some flexible function $h$ of a concept indicator and covariates accommodates many natural estimands of interest in the social sciences, though extensions to wider class of parameters that are associated with only asymptotically linear estimators will be discussed in detail in Section (ref).

examplecont{ex:rct} To infer average treatment effects on concept indicators, for treatment indicator $W_i$, fixed probability of treatment $\pi$ (bounded away from zero and one), and potential outcomes $Y_{ij}(w)$ under SUTVA, \[ H_{0,j}: E[Y_{ij}(1) - Y_{ij}(0)] = E[X_{ij}]=E\left[ \frac{W_{i} - \pi}{\pi(1-\pi)}Y_{ij}\right] = 0. \] In other words, we set $X_{ij} = h(W_i, Y_{ij}) := \frac{W_{i} - \pi}{\pi(1-\pi)}Y_{ij}$, a Horvitz–Thompson transformation for computing a difference in means.

We may form test statistics for such hypotheses as \[ T_{n,j} := \frac{1}{\sqrt{n}} \sum_{i=1}^n X_{ij}, \quad T_{n} := \frac{1}{\sqrt{n}} \sum_{i=1}^n X_{i}, \] and taking absolute values, as appropriate.

Though inference on all $p$ hypotheses communicates rich and interpretable information about our unstructured dataset, inference in this setting is also challenging because:

enumerate• The $p$ is large, so inference on all $j \in [p]$ parameters is inherently a multiple hypothesis testing problem. • Not only is $p$ large, but for most social science applications of interest $p \gg n$ (e.g., the experiment from Example (ref) contains a few hundred participants, but the concept embedding contains tens of thousands of concept indicators), so simultaneous inference on $Y_{ij}$ for each $j$ is a high-dimensional multiple hypothesis testing problem. • Each estimator, test statistic, or p-value formed from the $Y_{ij}$ for each $j \in [p]$ for the purposes of inference is plausibly statistically dependent on every other in a complicated way, ruling out multiple hypothesis testing approaches that assume independence or specific forms of dependence (e.g., PRDS). • Each $Y_i$ is sparse, and the intent of the analysis is discovery, so the desired form of selective (familywise) error control is ideally not too conservative.

What is needed, then, for principled simultaneous inference on parameters based on all $Y_{ij}$ is a dependence-agnostic, high-dimensional selective inference procedure with control over a generalized selective error rate. To be as well-powered as possible, we want to focus on testing procedures that employ resampling methods for estimating the true covariance across test statistics, such that conservative protection against worst-case dependence is not required. To keep things as general as possible, we also want to allow for only asymptotically valid test statistics and p-values, as is common in much of econometric analysis. In the following section, theory and corresponding statistical procedures are developed to achieve exactly these aims.

High-Dimensional Multiple Hypothesis Testing

We now state and prove the validity of statistical procedures that provide high-dimensional $k$-FWER control asymptotically, i.e., in the large sample limit, control of the probability of making $k$ or more rejections of true nulls under distribution $P$: \[ k\text{-FWER}_P := P\{\text{reject at least } k \text { hypotheses } H_{0,j}: j \in I(P)\} \] for $I(P)$ the set of indices of true null hypotheses under distribution $P$. Naturally, $k$-FWER control for $k=1$ is control of the FWER, and larger choices of $k$ power making more discoveries.

Background and Motivation

The fundamental challenge of high-dimensional statistical inference, characterized by settings where $p \gg n$, is that classic statistical theory no longer applies: classical asymptotic statistics treats $p$ as fixed in the relevant asymptotic thought experiment as $n$ grows large, whereas in high-dimensional settings the $p$ is growing (much) faster than $n$ in the relevant asymptotic thought experiment. As such, many classic results fail in high-dimensional settings, and new high-dimensional statistical theory must be developed---a literature that has grown popular in response to increasing interest in the statistical analysis of “big data.”\footnote{Concretely, challenges with developing statistical methods in high-dimensions include the lack of explicit limit distributions and complex dependence structures (see, e.g., chernozhukov_high-dimensional_2023 for more).} Multiple hypothesis testing (or selective inference, more broadly) in high-dimensions has received considerable attention in recent years, with methods now available that provide asymptotic high-dimensional control of the FWER under mild conditions belloni_high-dimensional_2018, and asymptotic high-dimensional control of the FDR under weak dependence across t-stats liu_phase_2014,belloni_high-dimensional_2018. Though valid in high-dimensions, these methods are quite conservative, or involve assumptions about the dependence structure of the data that are implausible in many settings (including the present one), respectively.

These existing high-dimensional multiple hypothesis testing results are supported by advances in high-dimensional central limit theory. Most existing high-dimensional CLT results naturally handle the case of controlling the FWER (which is facilitated by learning the quantiles of the distribution of the maximum coordinate of a scaled sum of independent random vectors in high-dimensions) but are invalid for studying the behavior of workhorse methods for controlling the $k$-FWER (which require modeling the distribution of the $k$-th largest coordinate of a scaled sum of independent random vectors in high-dimensions). Fortunately, recent work in high-dimensional central limit theory (i.e., ding_gaussian_2025) has developed new theoretical tools suited to understanding the behavior of exactly these objects in high-dimensions, and associated bootstrap procedures. We apply such tools to study and ultimately validate the asymptotic $k$-FWER control of the romano_control_2007 procedures in high-dimensions; we then extend these tools from ding_gaussian_2025, developing new high-dimensional CLTs that support more general high-dimensional multiple hypothesis testing procedures.

The following subsections state formal assumptions and results in high-dimensional central limit theory and multiple hypothesis testing. All formal results are stated as generally as possible, holding for data that is independent but not identically distributed (i.n.i.d.), which nests the setting of interest based on discovery from i.i.d. sampled unstructured data. Readers not interested in theoretical details and formal statements may skip to Section (ref).

Assumptions

Let $n \geq 3$ and $p \geq 3$. Let us further define \[ S_{n,j} := \frac{1}{\sqrt{n}} \sum_{i=1}^n (X_{ij} - E[X_{ij}]) , \quad S_n := \frac{1}{\sqrt{n}} \sum_{i=1}^n (X_i - E[X_i]). \] We will denote $\Sigma:=E\left[S_n S_n^\mathtt{T}\right]=n^{-1} \sum_{i=1}^n E[(X_i-E[X_i])(X_i-E[X_i])^\mathtt{T}]$, which simplifies to $\Sigma = E[(X_i-E[X_i])(X_i-E[X_i])^\mathtt{T}]$ under i.i.d. data.

Then we make the following assumptions.

assumption{M} Let $b_1>0$ and $b_2>0$ be some constants such that $b_1 \leq b_2$. For all $i\in [n],\, j \in [p]$, assume there exists sequence of constants $B_n>1$ (which may diverge) such that: (i) $E\left[\exp \left(\left|X_{i j}-E[X_{ij}]\right| / B_n\right)\right] \leq 2$; (ii) $b_1^2 \leq \frac{1}{n} \sum_{i=1}^n E[(X_{i j}-E[X_{ij}])^2]$; and (iii) $\frac{1}{n} \sum_{i=1}^n E[(X_{i j}-E[X_{ij}])^4] \leq B_n^2 b_2^2 $.

These are mild conditions on the tails and moments of the data $X_{ij}$, which are stated in a way that accommodates data that are i.n.i.d. Part (i) of Assumption (ref) simply requires that $X_{ij}$ be subexponential, or, equivalently, have an Orlicz norm in $\Psi_1$ bounded by $B_n$. The $B_n$ is indexed by $n$ to accommodate growing tail thickness within the subexponential regime as $p=p_n$ grows in the asymptotic thought experiment with i.n.i.d. data. Part (ii) of Assumption (ref) insists that the second (centered) moments of the data be bounded away from zero, appropriately stated for i.n.i.d. data. Part (iii) insists on bounded fourth (centered) moments in a similar fashion. Though these results are stated to handle i.n.i.d. data (in the interest of supporting applications beyond those considered in this paper), for i.i.d. data---the case considered in the applications of this framework---we may simply require for some fixed $B < \infty$ and fixed $b_1, b_2$ with $b_1 \leq b_2$: (i) $E\left[\exp \left(\left|X_{i j}-E[X_{ij}]\right| / B\right)\right] \leq 2$; (ii) $b_1^2 \leq E[(X_{i j}-E[X_{ij}])^2]$; and (iii) $ E[(X_{i j}-E[X_{ij}])^4] \leq B^2 b_2^2 $. Note that, in many use cases of interest based on concept indicators, the $X_{ij}$ will be both bounded (almost surely) and studentized, making these assumptions highly plausible.

remark[Degenerate features] Some fraction of concept indicators may be degenerate at zero for any given dataset (they are “dead” bricken2023monosemanticity). Assumption (ref) foreshadows that the high-dimensional CLTs underlying the upcoming high-dimensional multiple hypothesis testing results are only valid for understanding the $k$-th largest coordinate of the (in-population) non-degenerate coordinates of $S_{n}$, which we may label $S_{n}^+$. Of course, degenerate concepts in distribution $P$ are not of interest for the purposes of discovery, and so modeling the distribution of $S_{n,[k]}^+$ is precisely of scientific interest, in alignment with Assumption (ref). The researcher knows the in-sample degenerate coordinates, which is a superset of the in-population degenerate coordinates. We recommend dropping the in-sample degenerate coordinates to form a vector we may label $\check S_n$, and proceeding with analysis based on $\check S_n$ (noting that, by contrast, in finite samples, entries of $S_{n}^+$ may be zero, even if in the large sample limit they would not be). Whenever the $S_{n,j}^+$ for all non-degenerate $j$ are strictly positive with probability 1, then $\check S_{n,[k]} = S_{n,[k]}^+$ with probability 1 as well, and looking at the data to drop degenerate coordinates is both useful and benign. This commonly occurs: for example, when $|S_n^+|$ is of interest in two-sided testing and $X_{ij}$ has no atoms or $X_{ij}$ is binary (such as when $X_{ij}=Y_{ij}$).\footnote{So long as, technically, $nE[Y_{ij}]$ is not exactly an integer.} The Gaussian multiplier bootstrap quantity that will be discussed later on, under two-sided testing, will also satisfy this property under the data-conditional law. Otherwise, whenever $k$ or more coordinates of $S_{n}^+$ are strictly positive with very high probability, so too is $\check S_{n,[k]} = S_{n,[k]}^+$ with very high probability, and this will almost always be the case in conventional datasets of concept embeddings (that feature thousands of non-degenerate concepts).
assumption{R} Assume that $B_n^2 \log ^5(p n) = o(n)$.

This is the key rate condition needed for (uniform) Gaussian or bootstrap approximation error to go to zero asymptotically under the high-dimensional CLTs discussed in both the works of ding_gaussian_2025 and chernozhuokov_improved_2022. Rewritten, it says that \[ \frac{B_n^2 \log ^5(p n)}{n} = o(1) \] which permits $p$ growing very fast with $n$ in the asymptotic thought experiment. In fact, $p$ may be growing nearly exponentially in $n$, e.g., $p=e^{n^{1/6}}$ for fixed $B_n$. That this rate condition permits $p \gg n$ is an important positive result in high-dimensional central limit theory chernozhukov_central_2017,chernozhukov_high-dimensional_2023.

remark[State-of-the-art assumptions] Both Assumptions (ref) and (ref) are the state-of-the-art in the high-dimensional central theory from which this framework draws chernozhuokov_improved_2022,ding_gaussian_2025. Indeed, these conditions are the exact same as those used in both chernozhuokov_improved_2022 and ding_gaussian_2025. Notably, these assumptions do not require non-degeneracy of $\Sigma$, another important positive result in this literature chernozhukov_high-dimensional_2023. Under assumptions of non-degenerate $\Sigma$, improved rate conditions may be used in place of Assumption (ref).

High-dimensional $k$-FWER Control

We will focus on the Gaussian multiplier bootstrap as our high-dimensional bootstrap method, which is both computationally light and theoretically suited to the setting at hand. We define the Gaussian multiplier bootstrap quantity \[ S_n^B:=\frac{1}{\sqrt{n}} \sum_{i=1}^n \xi_i\left(X_i-\bar{X}_n\right), \] where $\xi_i \overset{iid}{\sim} N(0,1)$ and $\bar X_n := n^{-1} \sum_{i=1}^n X_i$. For the purposes of bootstrapping, we will denote the data-conditional probability measure $P^B(\cdot) := P(\cdot \mid X)$.

With these assumptions and notations in place, we can now state results that certify the application of the step-down algorithms of romano_control_2007 (Algorithms 2.1 and 2.2) for controlling the $k$-FWER in large samples in high-dimensions. Recall that Algorithm 2.1 of romano_control_2007 begins by assuming that the researcher can compute appropriate critical values given by $\hat{c}_{n, K}(1-\alpha, k)$ for any $K \subset \{1, ..., p\}$, which in the present case will be given by the $1-\alpha$ quantile of $S_{n,K,[k]}^B$ under $P^B$. The algorithm proceeds by first rejecting all $H_{0,j}$ such that $T_{n,j} > \hat{c}_{n, \{1,...,p\}}(1-\alpha, k)$, and ends if there are no rejections for any $j$, in which case all hypotheses are “accepted” (so far, this is equivalent to the entirety of a single-step procedure). If there are less than $k$ rejections, then the algorithm halts; otherwise, heuristically, the next step of the algorithm uses information about the previously rejected hypotheses to inform a choice of $K$ that yields a $\hat{c}_{n, K}(1-\alpha, k)$ that is then compared to the t-stats for the remaining hypotheses, leading to more rejections and another step (just like the previous) or halting if there are no more rejections. Algorithms 2.1 and 2.2 only differ on how they construct $K$ after the first step, with Algorithm 2.2 offering a more computationally streamlined approach. (For completeness, both algorithms are reproduced in full in Appendix Section (ref).)

restatable[High-dimensional $k$-FWER control for small $k$, one-sided]{theorem}{Smallk} Consider the method of Algorithm 2.1 in romano_control_2007 with test statistics $T_n$ of hypotheses $\{\tilde H_{0,j}\}_{j \in [p]}$ and critical values $\hat{c}_{n, K}(1-\alpha, k)$ given by the $1-\alpha$ quantile of $S_{n,K,[k]}^B$ under $P^B$. Assume that $k$ is fixed (i.e., not growing with $n,p$). Then under Assumptions (ref) and (ref), this method ensures: \begin{enumerate}[label=(\roman*)] • $\limsup_{n,p\to\infty} k\text{-FWER}_P \leq \alpha$. • If $\tilde H_{0,j}$ is false and $\theta_j(P) \gg B_n \sqrt{\log p /n}$, then the probability that the step-down method rejects $\tilde H_{0,j}$ tends to 1. • If $\min_{j \notin I(P)} \theta_j(P) \gg B_n \sqrt{\log p / n}$, then conclusions (i) and (ii) also hold when Algorithm 2.1 is replaced by Algorithm 2.2 of romano_control_2007. \end{enumerate}

This result is the natural high-dimensional analog of Theorem 3.1 (as well as Theorem 3.3) of romano_control_2007. Theorems 3.1 and 3.3 of romano_control_2007---as well as all formal results in that paper---are implicitly analyzed with fixed $p$ asymptotics. Theorem (ref) says that so long as $k$ is fixed (small in the relevant asymptotic sense), Algorithm 2.1 of romano_control_2007 may be used to asymptotically achieve control of the $k$-FWER at any desired level $\alpha$ even in high-dimensional settings, where $p$ can grow much faster than $n$, and moreover that Algorithm 2.1 is consistent against local alternatives that are not shrinking towards zero too quickly\footnote{Under Assumption (ref), note that $B_n \sqrt{\log p /n}=o(1)$.} (i.e., intuitively, in sufficiently large samples, we are able to appropriately reject a false null so long as $\theta_j(P)=E_P[X_{ij}]$ is not very close to zero). Result (i) of Theorem (ref) resembles result (i) of Theorem 3.1 of romano_control_2007, but where now $p$ is growing as fast as allowed by Assumption (ref); result (ii) of resembles result (ii) of Theorem 3.1 of romano_control_2007, but now $\theta_j(P) \gg B_n \sqrt{\log p /n}$ appears instead of $\theta_j(P) > 0$.\footnote{Formally, by $\theta_j(P) \gg B_n \sqrt{\log p /n}$ we mean that $\theta_j(P) /(B_n \sqrt{\log p /n}) \to\infty$.} Result (iii) says that the streamlined Algorithm 2.2 of romano_control_2007 also provides high-dimensional $k$-FWER control and is consistent against local alternatives when a natural (though perhaps not always granted) separation condition holds.

The proof of this result is available in Appendix Section (ref), and, as previously alluded to, relies on the insight that the low-dimensional $k$-FWER controlling algorithms of romano_control_2007 can be adapted to and analyzed in the high-dimensional setting using new results in high-dimensional central limit theory for order statistics chernozhuokov_improved_2022,ding_gaussian_2025. In the following theorem, we state the sibling procedure for high-dimensional control of the $k$-FWER in two-sided testing settings; it is the high-dimensional analog of Theorem 3.2 of romano_control_2007.

restatable[High-dimensional $k$-FWER control for small $k$, two-sided]{theorem}{Smallktwo} Consider the method of Algorithm 2.1 in romano_control_2007 with test statistics $|T_n|$ of hypotheses $\{ H_{0,j}\}_{j \in [p]}$ and critical values $\hat{c}_{n, K}(1-\alpha, k)$ given by the $1-\alpha$ quantile of $|S_{n,K}^B|_{[k]}$ under $P^B$. Assume that $k$ is fixed (i.e., not growing with $n,p$). Then under Assumptions (ref) and (ref), this method ensures: \begin{enumerate}[label=(\roman*)] • $\limsup_{n,p\to\infty} k\text{-FWER}_P \leq \alpha$. • If $ H_{0,j}$ is false and $|\theta_j(P)| \gg B_n \sqrt{\log p / n}$, then the probability that the step-down method rejects $H_{0,j}$ tends to 1. • If $\min_{j \notin I(P)} |\theta_j(P)| \gg B_n \sqrt{\log p / n}$, then conclusions (i) and (ii) also hold when Algorithm 2.1 is replaced by Algorithm 2.2 of romano_control_2007. \end{enumerate}
remark[Small $k$] To comment explicitly on why the procedures of Theorems (ref) and (ref) are only valid for small $k$, note that the sup-norm bootstrap or Gaussian approximation error in the high-dimensional CLTs for the $k$-th largest coordinate introduced in ding_gaussian_2025 only goes to zero if \[ \frac{k^8 B_n^2 \log ^5(p n)}{n} = o(1), \] meaning that $k$ must be fixed or must grow incredibly slowly: for all practical purposes, $k$ needs to be quite small. As such, the theoretical tools of ding_gaussian_2025 are unsuitable for making progress on FDP exceedance (FDX) control, as nontrivial FDX control at some $\gamma \in (0,1)$ (i.e., control of $P(\text{FDP} > \gamma)$) implies that $k$ is a non-negligible fraction of $p$ (i.e., $k$ grows nearly linearly with $p$).\footnote{More formally, as the proof of Theorem 4.1 in romano_control_2007 shows, $P(\text{FDP} > \gamma) = P(\text{false rejections} > k^*$), where $k^*$ is the smallest integer greater than $(\gamma /(1-\gamma))(p-|I(P)|)$.} A goal for future research would be to establish valid high-dimensional FDP exceedance control using further extensions of high-dimensional central limit theory. However, as will be observed in the empirical applications of this paper, small $k > 1$ are often sufficient to power meaningful discoveries in relevant social science settings at conventional choices of $\alpha$.
remark[Choice of $k$ and placebo tests] Though the high-dimensional central limit theory underlying Theorems (ref) and (ref) makes clear that $k$ must be “small,” it does not immediately suggest a particular (maximal) choice of $k$ for any finite sample. To make progress on choice of $k$, researchers may consider running placebo tests to help inform a reasonable choice of $k$. For example, in an experimental setting, researchers may generate placebo treatment assignments for units (independently of the data), and run the algorithm suggested by Theorem (ref) at increasing $ k \in \{1, 2, ...\}$. The first $k$ in this sequence that yields many rejections of the null relative to that $k$---which the researcher knows to be false rejections by placebo design---suggests that $k$ has become so large that control of the $k$-FWER has begun to fail. Such a procedure is used in the empirical example of Section (ref) to choose $k$.

Many Approximate Means and $k$-FWER Control

So far, we have validated the algorithms of romano_control_2007 in high-dimensions (with important caveats) using new theoretical tools for studying the $k$-th largest coordinate of $S_n$ from ding_gaussian_2025. However, these results were only based on hypotheses of means: the parameters of interest were written exactly as $\theta_j(P)=E_P[X_{ij}]$ for some interpretable $X_{ij}$, and estimators $\hat\theta_j$ were simply sample means (with $T_{n,j} = \sqrt{n}\hat\theta_j$). Though this is a common setting, researchers often want to test hypotheses for parameters which have estimators (and t-stats) that are not exactly expressible as averages of the observed data, yet “approximately” behave like means in large samples. Important examples include studentized t-stats,\footnote{Multiple hypothesis testing with studentized statistics is often important for the purposes of “balance”: that all hypothesis tests considered are similarly powered, contributing to selective error control similarly (see, e.g., romano_balanced_2010 for a more formal characterization of balance in this setting). For example, when $X_{ij} = Y_{ij} \in \{0,1\}$, more frequently occurring concept indicators are higher variance features, and we are typically interested in powering discoveries beyond just these frequently occurring concepts.} or when the parameters of interest are linear regression coefficients estimated with OLS. Thus, we also want $k$-FWER control for hypothesis testing based on “many approximate means” (MAM) in the parlance of belloni_high-dimensional_2018. To achieve this, we will extend the high-dimensional CLTs introduced in ding_gaussian_2025 (which are only valid for the “exact” means previously studied) to the setting of MAM, and thereby unlock high-dimensional $k$-FWER control for tests based on a wider class of estimators of parameters of practical interest.

To begin, note that many estimators $\hat\theta_j$ of parameters $\theta_j(P)$ will have the asymptotically linear representation \[ \sqrt{n}(\hat\theta_j - \theta_j(P)) = \frac{1}{\sqrt{n}}\sum_{i=1}^n \psi_{ij} + R_{n,j}, \quad \sqrt{n}(\hat\theta - \theta(P)) = \frac{1}{\sqrt{n}}\sum_{i=1}^n \psi_i + R_n, \] where $\hat\theta :=(\hat\theta_1, \hat\theta_2, ...,\hat\theta_p)^\mathtt{T}$, $\theta(P):=(\theta_1(P),\theta_2(P),...,\theta_p(P))^\mathtt{T}$, $\psi_i := (\psi_{i1},\psi_{i2},...,\psi_{ip})^\mathtt{T}$ is the collection of influence functions $\psi_{ij}$ at observation $i$ for $\theta_j(P)$, and the $R_{n,j}$ are small remainder terms (“linearization errors”). This setup nests the previous hypothesis testing setup based on “exact” means, for which $\theta_j(P) = E_P[X_{ij}]$ and \[ \mathbb{IF}(E_P[X_{ij}]) = X_{ij} - E_P[X_{ij}]=\psi_{ij}, \quad R_n = 0, \] where $\mathbb{IF}$ is the influence curve operator. Based on this observation, going forward we will generalize our earlier notation and simply write \[ S_n := \frac{1}{\sqrt{n}}\sum_{i=1}^n \psi_{ij}, \quad \hat S_n := \sqrt{n}(\hat\theta_j - \theta_j(P)) = S_n + R_n, \] letting $\Sigma:=n^{-1} \sum_{i=1}^n E[\psi_i \psi_i^\mathtt{T}]$, and treating Assumption (ref) as holding for $\psi_{ij}$.\footnote{Recall that, by definition, influence functions are mean zero, i.e., $E[\psi_{ij}]=0$.} With this new notation in place, we can now state a high-dimensional CLT that accommodates the case of small linearization errors $R_n$, the so-called setting of “many approximate means.”

restatable[High-dimensional CLT for the small $k$-max coordinate, MAM]{theorem}{Hdapprox} Let $\hat S_n := S_n + R_n$, and assume that $ \| R_n \|_\infty = o_P(1/\sqrt{\log p})$. Further assume that $k$ is fixed (i.e., does not grow with $n,p$). If Assumptions (ref) and (ref) hold, then as $n,p\to\infty$ \[ \sup _{t \in \mathbb{R}}\left|P\left( \hat S_{n,[k]} \leq t\right)-P\left(N(0,\Sigma)_{[k]} \leq t\right)\right| \to 0. \]

This theorem is both a natural generalization of Lemma A.6 of ding_gaussian_2025, which proves the special case for $R_n=0$, and an extension of Theorem 2.1 of belloni_high-dimensional_2018, which concerns uniform approximations over (hyper)rectangular events (like the max statistic falling below $t$) for “many approximate means.”\footnote{Note that the $k$-th largest coordinate falling below some $t$ for $k=1$ is an event where $\hat S_n$ falls in a (hyper)rectangular set, and for $k>1$ is a union of rectangles.} This theorem indicates that, for appropriately small $k$, in the large sample limit, the distribution of the $k$-th largest coordinate of a scaled sum of independent random vectors is (uniformly) well-approximated by the distribution of the $k$-th largest coordinate of a Gaussian random vector with variance matrix $\Sigma$ whether or not there is a small linearization error $R_n$. In Appendix Section (ref), we state an analogous theorem for the analogous (multiplier) bootstrap quantity under the relevant bootstrap law $P^B$.

Leveraging these high-dimensional CLTs for MAM, in Appendix Section (ref) we also state high-dimensional CLTs that explicitly allow for studentization, recognizing that studentization with estimated variances can be cast as analyzing a particular $\hat S_n$. These results can in turn be used to prove the following high-dimensional $k$-FWER result; a two-sided version can similarly be stated, but is omitted for brevity.

restatable[High-dimensional $k$-FWER control for small $k$, studentized, one-sided]{theorem}{Smallkstu} Define $\widehat{\Sigma}_{j j}:=n^{-1} \sum_{i=1}^n(X_{i j}-\bar{X}_{n, j})^2$ and $\widehat{\Lambda}:=\operatorname{diag}\{\widehat{\Sigma}_{11}, \ldots, \widehat{\Sigma}_{p p}\}$. Consider the method of Algorithm 2.1 in romano_control_2007 with test statistics $\hat\Lambda^{-1/2} T_n$ of hypotheses $\{\tilde H_{0,j}\}_{j \in [p]}$ and critical values $\hat{c}_{n, K}(1-\alpha, k)$ given by the $1-\alpha$ quantile of $(\hat\Lambda^{-1/2} S_{n}^B)_{K,[k]}$ under $P^B$. Assume that $k$ is fixed (i.e., not growing with $n,p$). Then under Assumptions (ref) and (ref), so long as $B_n=O(1)$, this method ensures: \begin{enumerate}[label=(\roman*)] • $\limsup_{n,p\to\infty} k\text{-FWER}_P \leq \alpha$. • If $\tilde H_{0,j}$ is false and $\theta_j(P) \gg \sqrt{\log p /n}$, then the probability that the step-down method rejects $\tilde H_{0,j}$ tends to 1. • If $\min_{j \notin I(P)} \theta_j(P) \gg \sqrt{\log p / n}$, then conclusions (i) and (ii) also hold when Algorithm 2.1 is replaced by Algorithm 2.2 of romano_control_2007. \end{enumerate}

Further note that, in the setting where the data is not just i.n.i.d. but i.i.d., as is true for the use case of interest based on concept indicators extracted from an unstructured dataset obtained through i.i.d. sampling, $B_n = O(1)$ is ensured (because it is constant).\footnote{It is also the case that we may safely only work on the event $\{\min _{j \leq p} \widehat{\Sigma}_{j j} \geq b_1^2 / 2\}$ in this theorem, so we incur no existence issues, as this event occurs with probability tending to 1 in the large sample limit given Assumption (ref).}

The evaluated influence functions $\psi_{ij}$ are typically not fully observed in the data, and instead must be estimated. This becomes important when bootstrapping, because the bootstrap must be based on quantities that we can actually compute with the observed data. The following theorem states that, when the influence functions are sufficiently well-estimated with observed data, the high-dimensional multiplier bootstrap remains effective.\footnote{Note that, even in the case of exact means, the influence function is not directly observed, because $E_P[X_{ij}]$ is not directly observed; nonetheless, the results proven earlier permit $E_P[X_{ij}]$ to be approximated with a sample average, as part of $S_n^B$. Theorem (ref) provides a higher level approximation condition to generalize this plug-in result.}

assumption{A} Assume that $\max_{j \in[p]} \left(n^{-1}\sum_{i=1}^n(\hat{\psi}_{i j}-\psi_{i j})^2\right)^{1/2}=o_P(1 / \log(pn)).$
restatable[High-dimensional bootstrap for the small $k$-max coordinate with estimated influence function]{theorem}{Hdbootinfluence} Define $$\tilde S_n^B := \frac{1}{\sqrt{n}}\sum_{i=1}^n \xi_i\hat\psi_i.$$ Further assume that $k$ is fixed (i.e., does not grow with $n,p$). If Assumptions (ref), (ref), and (ref) hold, then as $n,p\to\infty$ \[ \sup _{t \in \mathbb{R}}\left|P^B\left( \tilde S_{n,[k]}^B \leq t \right)-P\left(N(0,\Sigma)_{[k]} \leq t\right)\right| \xrightarrow[]{P} 0. \]

This theorem is the natural analog to Theorem 2.2 of belloni_high-dimensional_2018. It uses the same sufficient condition for the asymptotic approximation quality of the influence function---Assumption (ref)---as belloni_high-dimensional_2018 (equation 2.3), which is derived from an implication of the Borell-Sudakov-Tsirelson inequality.

The following result then leverages Theorem (ref), among others, to permit high-dimensional multiple hypothesis testing with $k$-FWER control in the MAM setting when influence functions must be estimated. Again, a two-sided result is omitted simply for brevity, and would introduce no new technical issues.

restatable[High-dimensional $k$-FWER control for small $k$, MAM, one-sided]{theorem}{Smallkaif} Consider the method of Algorithm 2.1 in romano_control_2007 with test statistics $ \sqrt{n}\hat\theta_j$ of null hypotheses $\tilde H_{0,j}: \theta_j(P) \leq 0$ for all $j \in [p]$, and critical values $\hat{c}_{n, K}(1-\alpha, k)$ given by the $1-\alpha$ quantile of $\tilde S_{n,K,[k]}^B$ under $P^B$. Assume that $k$ is fixed (i.e., not growing with $n,p$) and that we may write for all $j \in [p]$ \[\sqrt{n}(\hat\theta_j - \theta_j(P)) = \frac{1}{\sqrt{n}}\sum_{i=1}^n \psi_{ij} + R_{n,j}\] where $ \sup_{j\in [p]} |R_{n,j}| = o_P(1/\sqrt{\log p})$ and $\psi_{ij}$ is the influence function for $\theta_j(P)$ evaluated at observation $i$. Then under Assumptions (ref), (ref), and (ref), this method ensures: \begin{enumerate}[label=(\roman*)] • $\limsup_{n,p\to\infty} k\text{-FWER}_P \leq \alpha$. • If $\tilde H_{0,j}$ is false and $\theta_j(P) \gg B_n\sqrt{\log p /n}$, then the probability that the step-down method rejects $\tilde H_{0,j}$ tends to 1. • If $\min_{j \notin I(P)} \theta_j(P) \gg B_n \sqrt{\log p / n}$, then conclusions (i) and (ii) also hold when Algorithm 2.1 is replaced by Algorithm 2.2 of romano_control_2007. \end{enumerate}

Theorem (ref) embeds a practical bootstrapping procedure for estimators of parameters that are only asymptotically linear: derive the influence functions for the parameters of interest, estimate them if they are not directly observed, and then implement a multiplier bootstrap based on the data-conditional law of $\tilde S_n^B$.

An important application of these results is the empirical setting where each $\theta_j(P)$ is a linear regression coefficient based on a separate regression of a concept indicator on some regressor(s) of interest, estimated with a corresponding OLS coefficient $\hat\theta_j$. Here, $\sqrt{n}(\hat\theta_j -\theta_j(P))$ will have only an asymptotically linear representation, and the data-evaluated influence functions must be estimated. In the following example, we walk through the derivation of the bootstrapping procedure for testing many regression coefficients.

example[Bootstrap for many linear regression coefficients] Consider the dataset given by $\{(Y_i,T_i,D_i)\}_{i=1}^n$, drawn i.i.d. from distribution $P$, where $Y_i=(Y_{i1},\ldots,Y_{ip})^{\mathtt T}\in\mathbb R^p$ are concept indicators, $T_i\in\mathbb R$ is a treatment of interest, and $D_i\in\mathbb R^d$ are controls for some small (fixed) $d$ that includes a constant. The researcher is interested in making discoveries based on the coefficient on $T$ from the linear regression of $Y_j$ on $(T, D)$, which we denote $\theta_j(P)$, and which is formally defined via the optimization problem \[ (\theta_j(P),\vartheta_j(P)) := \arg\min_{b,g} E_P[(Y_j - T b - D^\mathtt{T} g)^2]. \] By the Frisch-Waugh-Lovell theorem, we can write \[ H_{0,j}: \theta_j(P) = \Omega(P)^{-1} E_P[\widetilde T\,\widetilde Y_j]=0 \] where $L_P[\cdot\mid\cdot]$ denotes the population linear projection operator, and by defining $\widetilde T := T - L_P[T\mid D]$ and $\widetilde Y_j := Y_j - L_P[Y_j\mid D]$, letting $\Omega(P) := E_P[\widetilde T^2]$ for shorthand. Defining the population residual $\widetilde U_j := \widetilde Y_j - \theta_j(P) \widetilde T$, the influence function for $\theta_j(P)$ at observation $i$ is given by \[ \psi_{ij} := \Omega(P)^{-1}\widetilde T_i \widetilde U_{ij}. \] For the expectation under the empirical distribution $E_n[\cdot]$ (i.e., the sample average), we can define feasible sample analogs as \[ \begin{gathered} \widehat{\gamma}_T:=E_n[D_i D_i^{\mathtt{T}}]^{-1}E_n[ D_i T_i],\qquad \widehat{\gamma}_{Y_j}:=E_n[D_i D_i^{\mathtt{T}}]^{-1}E_n[D_i Y_{i j}],\\ \widehat{\widetilde{T}}_i:=T_i-D_i^{\mathtt{T}} \widehat{\gamma}_T, \qquad \widehat{\widetilde{Y}}_{ij}:=Y_{ij}-D_i^{\mathtt{T}} \widehat{\gamma}_{Y_j},\\ \widehat\Omega:=E_n[\widehat{\widetilde T}_i^2], \qquad \widehat\theta_j:=\widehat\Omega^{-1}E_n[\widehat{\widetilde T}_i\,\widehat{\widetilde Y}_{ij}],\qquad \widehat{\widetilde U}_{ij}:=\widehat{\widetilde Y}_{ij}-\widehat\theta_j\widehat{\widetilde T}_i, \qquad \widehat\psi_{ij}:=\widehat\Omega^{-1}\widehat{\widetilde T}_i \widehat{\widetilde U}_{ij}, \end{gathered} \] noting that $\hat\theta_j$ is just the usual OLS coefficient. We may now write out the asymptotically linear representation of each of the $j\in [p]$ regression coefficients as \[ \sqrt{n}(\hat\theta_j - \theta_j(P)) = \frac{1}{\sqrt{n}}\sum_{i=1}^n \psi_{ij} + R_{n,j}. \] The appropriate law for facilitating $k$-FWER control via Theorem (ref) is then given by \[ \tilde S_n^B := \frac{1}{\sqrt{n}}\sum_{i=1}^n \xi_i\hat\psi_i = \widehat\Omega^{-1}\frac{1}{\sqrt{n}}\sum_{i=1}^n \xi_i \widehat{\widetilde T}_i \widehat{\widetilde U}_{i} \] under $P^B$, where $\widehat{\widetilde U}_{i}$ stacks all $\widehat{\widetilde U}_{ij}$ as a vector. That $\hat\psi_i$ is a sufficiently good approximation of $\psi_i$ and that $R_n$ is sufficiently small (in large samples) is proven in Propositions (ref) and (ref), as stated in Appendix Section (ref), under certain sufficient conditions (e.g., all data is subgaussian tailed).
commentIn some instances, standard deviations computed in sample may be zero, because of the sparsity with which the dictionary features activate. As such, we may consider a modified estimator $\tilde \Lambda$ where \[ \tilde \Lambda_{jj} := \max\left\{ \hat \Lambda_{jj}, \underline{\lambda} \right\} \] for some $\underline{\lambda} > 0$. We may choose $\underline{\lambda}$ based on a small minimum detectable effect of interest, with the only consequence being that tests of smaller effects than this minimum are plausibly inconsistent. For example, for $X_{ij} = Y_{ij} \in \{0,1\}$, consider setting the minimum detectable probability of dictionary feature activation of interest as $\tau > 0$; then one may reasonably set $\underline{\lambda} = \tau (1- \tau)$.

Generating and Evaluating Automatic Interpretations

By applying the statistical procedures discussed in Section (ref), the researcher obtains a selected set of rejected null hypotheses, or “discoveries.” Concretely, these discoveries are recorded as set of coordinates $\hat J \subseteq [p]$ indexing the rejected nulls.

Though the machine learning theory that motivates sparse dictionary learning supposes that each hypothesis $\hat \jmath \in \hat J$ has single, human-interpretable meaning, dictionary learning methods themselves do not generate natural language descriptions of what each human-interpretable meaning may be. We now discuss a set of procedures for both automatically generating and evaluating the quality of such natural language feature descriptions from a statistical inference perspective.

remark[Why do we need automatic interpretation?] Human judgment alone could be deployed as a means of giving natural language descriptions to concept indicators, e.g., the researcher could inspect many of the text examples that activate (highly) for a given SAE feature in order to arrive at a plausible English language description for the concept signaled by that feature. However, such a process suffers from (1) problems related to motivated reasoning and data snooping (e.g., the researcher could simply come up with a concept description that is reverse engineered from a conclusion they wanted to find in the data); and (2) problems related to the scale of interpretation and associated interpretation effort/cost (e.g., there may be many hundreds of concepts discovered, all requiring interpretation; interpreting all $p$ concepts would require even more labor). Automatic interpretation pipelines solve both of these problems. They automate the process of creating and scoring natural language descriptions of concept indicators by leveraging a dedicated “explainer LLM” to reason about interpretations cheaply at scale.\footnote{The terminology “explainer LLM” is chosen to be consistent with related works in the AI interpretability literature, e.g., jiang_interpretable_2025.} (Of course, the researcher can look directly at the highly activating unstructured data points to develop their own interpretations, simultaneously.)

Automatic Generation

We adopt the autointerp generation strategy pioneered in bricken2023monosemanticity for this framework: we pass text samples from the researcher's unstructured dataset that activate highly on a particular feature to an explainer LLM---annotating these exemplar texts with information about where and how strongly the feature activates---and ask for an English description of the possible meaning of the feature.\footnote{This exact generation strategy, which only looks at “top activating” examples, is also used in and argued for in, e.g., choi2024automatic.}

By drawing on text samples from the researcher's unstructured dataset, this generation strategy is tailored to the distribution $P$. Though off-the-shelf feature descriptions are sometimes available for pretrained dictionary learning models, these descriptions were necessarily learned on a different distribution of data (e.g., some large and opaque swath of the internet). As such, these descriptions are liable to be of lower quality than those based on the researcher's own dataset. The specific prompts implementing this strategy for the empirical applications considered later are detailed in Appendix Section (ref).

examplecont{ex:rct} The null hypothesis of no average treatment effect is rejected for concept 4. Using the raw SAE feature activations (necessary to compute in the process of evaluating the Dict function), the researcher is able to rank text responses from experimental participants by largest activations for concept 4 on any token (“max pooling”). The researcher takes the top 10 highest activating text examples and annotates each of them with brackets to delimit the highest activating token(s) in each response. These exemplar responses are then passed to an explainer LLM that reasons about and then outputs the concise English phrase “discussion of politics” to describe concept 4.

Automatic Evaluation: A Statistical Inference Perspective

Autointerp pipelines ideally not only specify how to generate descriptions for concept indicators, but also evaluate the quality of these descriptions. We now consider a new statistical formalization of a mainstream autointerp evaluation technique called “detection scoring.”

The fundamental principle of detection scoring is that “an interpretation should serve as a binary classifier distinguishing activating from nonactivating contexts” paulo_automatically_2024. In other words: description of a concept can be treated like a classification problem, one where the natural language description itself (harnessed in a LLM prompt) is the classifier, the input is a specific text example, and the label to be predicted is derived from the ground truth of the presence or absence of a SAE feature's activation in that text example (which is always known to the researcher using their Dict function). Under this framing, binary classifier performance metrics like accuracy, precision, and recall characterize the quality of the description.

To formalize this autointerp strategy, consider a “held-out” unstructured dataset of size $m$ drawn independently from the same distribution of interest $P$. (This held-out dataset may be created, for example, by partitioning an original unstructured dataset drawn i.i.d. from $P$ into two disjoint datasets via indices chosen at random, independently of all the data.) Let us denote the set of indices of unstructured data points in the main sample as $\mathcal{I}^\textrm{estim}$, and those in the held-out dataset as $\mathcal{I}^\textrm{eval}$.\footnote{In empirical settings featuring a randomized experiment, for sake of preserving power, it may be useful for researchers to create the held-out dataset split from the larger treatment arm, if one exists.}

Using the data points with indices in $\mathcal{I}^\textrm{estim}$, we extract exemplar unstructured data points to be used to create an autointerp explanation of concept $j$, which we denote $\hat\eta_j$. For all $i \in \mathcal{I}^\textrm{eval}$, we define the random variable $\hat Y_{ij} \in \{0,1\}$ as \[ \hat Y_{ij} := \texttt{CLS}(Z_i, \hat \eta_j), \] where the function $\texttt{CLS}(\cdot)$ is a fixed choice of LLM and classification prompt fed both an unstructured data point $Z_i$ and the learned autointerp explanation $\hat\eta_j$. The $\hat Y_{ij}$ is the classifier output, a prediction of whether or not the unstructured data point $Z_i$ contains a concept with description $\hat \eta_j$.\footnote{Note that, even if a sample splitting implementation of this autointerp evaluation strategy is pursued, cross-fitting is inappropriate, because there is little reason to think, given the sample sizes available to researchers, that $\hat\eta_j$ would be stable across folds. As such, when sample splitting is employed, an uneven split is recommended, e.g., only 10% of the data is used for the held-out evaluation sample.}

For the purposes of evaluating autointerp descriptions, we define the “accuracy score” (or “A-score”) estimand as \[ \theta_j^\textrm{acc}(\hat\eta_j;P) := E_P\left[S_{ij} \mid \hat\eta_j \right] = P(\{Y_{ij} = \hat Y_{ij}\} \mid \hat\eta_j) = P\left(\{Y_{ij} = \texttt{CLS}(Z_i, \hat\eta_j)\} \mid \hat\eta_j\right) \] where $S_{ij} := \mathbf{1}\{Y_{ij} = \hat Y_{ij}\}$. Conditioning on $\hat\eta_j$ in the estimand of interest should be thought of as evaluating a fixed autointerp description. Importantly, as the notation emphasizes, the $\theta_j^\textrm{acc}(\hat\eta_j; P)$ tells us exactly the quality of the description $\hat\eta_j$ for the population of interest $P$, unlike other autointerp evaluation methods for which the underlying distribution that implicitly defines the evaluation score is highly opaque (c.f., simulation scoring bills2023language).

We can easily construct a plug-in estimator of the autointerp accuracy (A-score) estimand; its unbiasedness, consistency, and asymptotic normality under the $\hat\eta_j$-conditional law are established in the following proposition.

restatable[Conditional inference on A-score]{proposition}{Ascoreest} Let $m:=|\mathcal{I}^\textrm{eval}|$. Define the A-score estimator for concept indicator $j$ as \[ \hat \theta_j^\textrm{acc} := \frac{1}{m} \sum_{i \in \mathcal{I}^\textrm{eval}} S_{ij} ,\quad S_{ij}= \mathbf{1}\{Y_{ij} = \hat Y_{ij}\} = \mathbf{1}\{Y_{ij} = \operatorname{\texttt{CLS}}(Z_i, \hat \eta_j)\}. \] Assume that $P(S_{ij}=1 \mid \hat \eta_j) \in [\delta, 1-\delta]$ for some $\delta$ bounded away from 0 and 1, almost surely. Then we have that $E_P\left[\hat \theta_j^\textrm{acc} \mid \hat\eta_j \right]=\theta_j^\textrm{acc}(\hat\eta_j;P)$ almost surely, and that as $m \to \infty$ \[ \sup_{t \in \mathbb{R}} \left| P\left(\sqrt{m}( \hat \theta_j^\textrm{acc} - \theta_j^\textrm{acc}(\hat\eta_j;P)) \leq t \mid \hat\eta_j \right) - P(N(0, \operatorname{Var}(S_{ij} \mid \hat\eta_j)) \leq t \mid \hat\eta_j) \right| \xrightarrow[]{a.s.} 0. \]

It is immediate from Proposition (ref) then that the empirical variance estimator $\hat V_j^\text{acc}:= m^{-1}\sum_{i \in \mathcal{I}^\text{eval}} (S_{ij}-m^{-1}\sum_i S_{i j} )^2$ permits construction of valid large sample frequentist confidence intervals for $\hat\theta_j^\text{acc}$ under the $\hat\eta_j$-conditional law of the data.

examplecont{ex:rct} The null hypothesis of no average treatment effect is rejected for concept 4, and the autointerp generation strategy implemented by the researcher describes concept 4 as activating on “discussion of politics.” To evaluate the quality of this description, the researcher runs a LLM-based classifier over a held-out dataset of experimental participant text responses. The LLM classifier is given the prompt “Does the following text example indicate discussion of politics? Return 1 or 0 for yes or no:” followed by a text response from the held-out dataset. The researcher then uses the output of the LLM classifier to compute the A-score, which measures the accuracy of this autointerp description-based classifier at predicting the presence of concept 4 in the held-out dataset.
remark[Oracle detection scores] Towards better understanding the information communicated to the researcher by the A-score, it is useful to consider the function $\texttt{CLS}^*(\cdot)$, which we can define implicitly as \[ \bar\theta_j^\textrm{acc}(\hat\eta_j;P):=\sup_{f \in \mathcal{F}} P\left(\{Y_{ij} = f(Z_i, \hat\eta_j)\} \mid \hat\eta_j\right) = P\left(\{Y_{ij} = \texttt{CLS}^*(Z_i, \hat\eta_j)\} \mid \hat\eta_j\right), \] where $\mathcal{F}$ is the space of all viable classification prompts and LLMs able to execute them. That is, $\texttt{CLS}^*(\cdot)$ is the “oracle” detection scoring model setup, which when evaluated produces the highest interpretative accuracy score possible under a given description $\hat\eta_j$ under distribution $P$. Arguably, $\bar\theta_j^\textrm{acc}(\hat\eta_j;P)$ is more informative to the researcher than $\theta_j^\textrm{acc}(\hat\eta_j;P)$ about the quality of $\hat\eta_j$, as $\bar\theta_j^\textrm{acc}(\hat\eta_j;P)$ purely quantifies the quality of $\hat\eta_j$, without the wedge introduced by classification failures stemming from $\texttt{CLS}(\cdot)$ itself. Because any feasible choice of $\texttt{CLS}(\cdot)$ yields \[ P\left(\{Y_{ij} = \texttt{CLS}(Z_i, \hat\eta_j)\} \mid \hat\eta_j\right) \leq P\left(\{Y_{ij} = \texttt{CLS}^*(Z_i, \hat\eta_j)\} \mid \hat\eta_j\right) \] one can reasonably interpret $\hat \theta_j^\textrm{acc}$ as a conservative (under)estimate of the more informative quantity $\bar\theta_j^\textrm{acc}(\hat\eta_j;P)$.

We may also define estimands for autointerp precision and recall (the “P-score” and “R-score,” respectively):

align*[align* omitted — 455 chars of source]

As above, these estimands may be consistently estimated with appropriate asymptotically valid confidence intervals under the $\hat\eta_j$ conditional law using the natural plug-in estimators and the delta method. That said, for the purposes of the empirical applications in the following section, we focus on estimating the A-score, following common practice in the relevant literature paulo_automatically_2024.

Further, recall that we are not typically interested in autointerp evaluation for a fixed $j$, but instead for a particular data-driven choice $\hat\jmath \in \hat J$ suggested by a high-dimensional selective inference procedure.\footnote{Autointerpretation typically incurs non-negligible compute costs, and, as such, researchers will often only be interested in forming meaningful autointerpretations for the subset of discovered concepts.} However, so long as discovery is also performed using the main estimation sample indexed by $\mathcal{I}^\textrm{estim}$, the researcher may treat each $\hat \jmath$ as fixed for the purposes of evaluation inference. Formally, for all $j \in [p]$,

align*[align* omitted — 200 chars of source]

where the last equality follows because, under the proposed independently held-out data scheme, \[ [\mathbf{1}\{Y_{ij} = \texttt{CLS}(Z_i, \hat\eta_j)\} \perp\!\!\!\perp \hat \jmath] \mid \hat\eta_j. \]

remark[Inference on the best interpretation] If a researcher is interested in inference on the best performing autointerp description across multiple generation prompts or strategies---say among a set $\{\hat\eta_{j,1}, \hat\eta_{j,2}, \dots\}$ for any $j$ of interest---then the evaluation performance should be estimated with appropriate “inference on winners” corrections such as those discussed in andrews_inference_2024.
remark[Evaluation of human-generated descriptions] It is worth noting that this autointerp evaluation strategy can also serve as a general-purpose description evaluation technology, e.g., for carefully crafted human descriptions of concept indicators in settings where there is a small number of discoveries. The scientific benefit of evaluation in such a setting is also clear: even if there is concern that human-made descriptions are reverse-engineered to support a presupposed conclusion, any description (human made or machine made, motivated or not) that is not supported by the data will be surfaced through the evaluation scores. In this way, description evaluation can be seen as even more critical when automatic methods are not employed by the researcher to generate descriptions.

Pseudocode and Code

In Appendix Section (ref), we state the entire framework for interpretable discovery from unstructured data as a single pseudocode algorithm. We focus on an implementation with text data, sample splitting, two-sided testing, single-step hypothesis testing (i.e., just the first step of Algorithm 2.1 of romano_control_2007), and autointerp evaluation via accuracy scoring (A-scores). Full implementation code is also available: see the link in Appendix Section (ref).

Empirical Applications

As an illustration of the framework proposed in this paper, we now reanalyze two recent works in empirical economics that pursue discovery from unstructured data, and show how new, principled, and interpretable discoveries may be made using the same exact data at low cost.

Reanalyzing bursztyn_justifying_2023

Mirroring Example (ref), bursztyn_justifying_2023 study how the provision of “social cover” affects willingness to publicly dissent to socially stigmatized causes, and the perception of this dissent. As a key application of their theory, bursztyn_justifying_2023 run an information treatment experiment (Experiment 2) in which participants are told they have been matched with another (fictional) respondent that chose to join a campaign to oppose defunding the police (a plausible expression of dissent in liberal American politics at the time the experiment was conducted), and then show the participant a tweet that the matched respondent is said to have agreed to publicly post. This tweet has been randomized to either provide or not provide some “social cover” for the fictional matched respondent's plausibly contentious statement, namely whether or not the tweet indicated that the matched respondent joined the campaign before reading an evidence-based article in support of it (the no social cover condition) or after (the social cover condition).\footnote{See Online Appendix Figure B.3 of bursztyn_justifying_2023 for a schematic of the experimental design.}

An key outcome collected by this experiment is the participant's open-ended text response to the question “Why do you think your matched respondent chose to join the campaign to oppose defunding the police?” This open-ended text response is meant to capture the causal effect of social cover on the perception of dissent. As such, this is a setting where discovery is of (at least partial) interest: ideally, to form a holistic understanding of how social cover affects perception of dissent, we do not want to ex-ante pre-specify what aspects of perception we are interested in detecting, and want to instead discover any interpretable systematic differences that exist in the text responses across treatment and control groups. In order to make progress on this analysis of unstructured data, bursztyn_justifying_2023 compute a Pearson's $\chi^2$ statistic for all phrases of up to three words per gentzkow_what_2010, which they use as an index to rank the phrases that are most differentially expressed in each condition's open-ended responses.\footnote{Specifically, it is a Pearson's $\chi^2$ statistic for a null hypothesis that the propensity to use a given phrase is equal across conditions, per gentzkow_what_2010. We interpret these statistics as simply forming an index, however, because the results of these hypothesis tests are not reported, and no multiple hypothesis testing corrections are implemented.} The interpretability of the results from this analysis, is, naturally, hindered by the coarseness of the featurization of the text as three word phrases, as well as the fact that there are no obvious estimands or inferential guarantees. The only qualitative conclusion gleaned by bursztyn_justifying_2023 from this quantitative exercise is: “we find that respondents in the Cover condition are more likely to use phrases related to the article or the associated evidence—for example, `article,' `read,' `convincing,' or `increase in crime.'\,”\footnote{The Online Appendix Table B.11 of bursztyn_justifying_2023 contains the top ten characteristic phrases in each condition based on the $\chi^2$ index.}

Can we make more (and more interpretable) discoveries with the framework proposed in this paper? To investigate, we implement a $\texttt{Dict}$ function based on the Gemma Scope autoencoders for Google's Gemma 2 2B LLM lieberum_gemma_2024, and specifically the SAE trained on the residual stream of layer 12.\footnote{It is a well-documented observation that SAEs trained on the middle layer of a LLM seem to capture “coarser” concepts and features (see discussion in, e.g., templeton2024scaling), which is plausibly more suitable for analyzing broad themes in text.} We implement an autointerpretation pipeline as discussed in Section (ref) using as an explainer LLM the lightweight, open-source reasoning model GPT-OSS 20B openai2025gptoss120bgptoss20bmodel. A held-out evaluation dataset is created by sample splitting, using just 10% of the original dataset in bursztyn_justifying_2023.

The hypotheses of interest are, for all $j \in [p]$, \[ H_{0,j}: E[Y_{ij}(1) - Y_{ij}(0)]= E\left[\frac{W_{i} - \pi}{\pi(1-\pi)}Y_{ij}\right] := E[X_{ij}] = 0, \] for $\pi = 0.5$, and where the event $\{W_i = 1\}$ indicates the social cover condition. This yields a dataset $\{X_{ij}\}_{i\in[n], j\in[p]}$, where $n = 930$ and $p = 11901$ (total non-degenerate features---see Remark (ref)). We implement the simple one-step version of the statistical procedure in Theorem (ref) as per Algorithm (ref), which attains a lower bound on the power achieved by the step-down procedure.

We may first consider whether or not, as anticipated, simply controlling the FWER would be too conservative for the purposes of discovery in this setting. To explore this, we set $k=1$ and $\alpha=0.05$ (equivalently, we implement the high-dimensional FWER controlling procedure of belloni_high-dimensional_2018), and we obtain only two significant discoveries: a 13.2 percentage point causal effect on the presence of concept 3518 and a 13.4 percentage point causal effect on the presence of concept 3426. The autointerprations of these concepts paint a picture that corroborates the original findings of bursztyn_justifying_2023: concept 3426 is described as activating on mentions of “defunding the police associated with increased violent crime” and concept 3518 on “mention of an article as a source of evidence.” However, we only discover two concepts out of approximately 12,000, and not much new insight relative to the original analysis---even if the discovery and interpretation of these concepts was transparent, principled, and automatic.\footnote{The off-the-shelf Gemma Scope autointerpretations echo the same: concept 3426 is described as activating on “references to news articles and reports” and concept 3518 on “references to articles.”}

figure[figure omitted — 171 chars of source]

In Figure (ref), we make discoveries now setting $k=5$, still with $\alpha = 0.05$.\footnote{See Remark (ref) for how this choice of $k$ was made.} As can be seen, even with control of $k$-FWER at $k=5$, we obtain 9 rejections of the null, or 9 discoveries. By large sample guarantee of Theorem (ref), the probability of making 5 or more false rejections with this procedure is at most 5%, and so, informally speaking, we should be confident that at least 5 of these 9 discoveries are genuine.\footnote{To get a better sense for which discoveries are plausibly true, one might consider re-running the algorithm of Theorem (ref) at increasing values of $k$, and noting at which $k$ each hypothesis enters the selected set.} The previous concepts discovered with (1-)FWER control appear, as well as others absent from the (1-)FWER controlled analysis and the original analysis in bursztyn_justifying_2023. Though there are redundant concepts discovered, in totality these discoveries contribute to a few entirely new characterizations of the causal effect of social cover on dissent. For example, concept 1622 responds to “Prestigious author source reference” and increases by nearly 8 percentage points in the treatment group, indicating that---exactly as bursztyn_justifying_2023 posit, and run an entirely separate experiment to further probe---the social cover mechanism plausibly relies on the perceived credibility of the article cited. Concept 1355 responds to “modal verbs expressing uncertainty or possibility,” and its probability of activation is decreased by nearly 12 percentage points in the treatment group. This supports the idea that the causal effect of not having social cover is, to first order, generating uncertainty and speculation about dissent, as opposed to an alternative scenario where participants coalesce around a single narrative in the absence of cover. This insight is also congruous with the observation that few large negatively signed causal effects are discovered. The A-scores for all discoveries range from 0.98 to 0.66, indicating that some descriptions of discoveries deserve more skepticism than others, and that there is much scope for improving autointerp generation beyond small, open-source reasoning LLMs with the prompting strategy of Appendix Section (ref).

Importantly, all of this analysis was performed automatically, in approximately one hour, on a single A100 GPU in Google Colab (for free, as a student or university affiliate). Though this analysis was not preregistered, it is possible to cheaply and automatically replicate these results, and assess their sensitivity to different autointerpretation strategies, concept embedding implementations, and testing hyperparameters. Notably, the space of possible concepts was not pruned in any way prior to this analysis, nor were they weighted in some way towards topics related to articles; the analysis was conducted without human intervention, with no room for data snooping.

Reanalyzing stantcheva_why_2024

Mirroring Example (ref), stantcheva_why_2024 conducts surveys on representative samples of the United States population investigating attitudes towards inflation. These surveys include open-ended text responses to a variety of questions, with the goal of discovering attitudes and opinions that the researcher should not (for concern about priming respondents) or could not (for lack of imagination) manually pre-specify.\footnote{See, e.g., haaland_understanding_2024 for more discussion of the benefits of open-ended survey questions for understanding economic behavior.}

An important open-ended prompt that stantcheva_why_2024 solicit an answer to is “High inflation is caused by...” In the original analysis, stantcheva_why_2024 codes open-ended responses with a discretionary keyword-based method introduced in ferrario_eliciting_2022 (which suggests the researcher create a list of topics and associated keywords that denote their presence via anything from “manual to semi-supervised or unsupervised” means). The result of applying this analytic technique to the open-ended responses to the above prompt yields Figure 3 in stantcheva_why_2024, which finds, among other things, that mentions of “Biden and the administration,” “Greed,” “Monetary policy,” “Fiscal policy,” “War and foreign policy,” “Demand vs supply,” and “Supply-side mechanisms (other than input prices)” appear in more than 5% of all responses.

Do we make the same discoveries automatically if we apply the proposed framework instead? Do we make more and new interpretable discoveries? To find out, we use the same $\texttt{Dict}$ function and explainer LLM as in Section (ref), but reduce the space of SAE features by half, filtering out all features that, in the corpus on which the Gemma Scope SAEs were trained, had greater than median empirical token activation frequency. This choice of dimensionality reduction is meant to screen out features that activate in many texts across domains, e.g., features related to grammatical aspects of text, which was handled implicitly in the previous causal analysis by virtue of differencing. (Ideally, such a dimensionality reduction choice would be preregistered, to prevent cherry-picking by filtering, but could also easily be analyzed for sensitivity.) Again we split the sample randomly to form a held-out evaluation dataset consisting of 10% of the original dataset analyzed in stantcheva_why_2024.

We now test hypotheses $H_{0,j}: E[X_{ij}] := E[Y_{ij}] = 0$, i.e., we test the probability of a given concept being discussed in the population of open-ended responses to the above prompt in the United States. (Note that, in this setting, rejecting $E[Y_{ij}] = 0$ is less of interest than the dual of this testing procedure, which generates point estimates for the probability of discussion along with generalized simultaneous confidence intervals.) This yields a dataset $\{X_{ij}\}_{i \in [n], j\in [p]}$ with $n=453$ and $p=3826$ non-degenerate features---still a very high-dimensional inference problem. We now apply 5-FWER control with $\alpha = 0.05$ based on the result in Theorem (ref) per Algorithm (ref), and we make 744 interpretable discoveries.

In Figure (ref), we report point estimates and simultaneous confidence intervals for all discoveries with a generalized simultaneous lower confidence bound above 20%. Notably, we discover, automatically, with no manual pruning or interpretation required, many sensible commonplace themes related to the topic at hand: concept 9804 responds to discussion of the economy; concept 13447 responds to mentions of inflation; and concept 13574 activates on broad discussion of negative economic conditions. However, we also recover many of the topics from the original analysis in stantcheva_why_2024 without any human discretion or intervention: concept 4192 responds to discussion of monetary policy; concept 14747 responds to discussion of supply-side issues; concept 104 activates on discussion of excess and overabundance; concept 8316 responds to fiscal policy mentions; and concept 5719 responds to economic criticisms leveled at the current government. Undoubtedly, many other topics are recovered as we look deeper down the list of discoveries, beyond these 11 largest. However, even in this top 11, we learn that there is a great deal of pure uncertainty expressed in these open-ended responses based on the activation of concept 4794, and also that macroeconomic indicators get frequently discussed by virtue of concept 11036.

figure[figure omitted — 193 chars of source]

Once again, this procedure was inexpensive and fast to implement, and could be quickly replicated and analyzed for sensitivity by any other researcher. Beyond filtering out greater than median activating concepts, no other choices were made to delimit the space of possible discoveries, or focus them in some way on the space of economic- or inflation-relevant concepts; the exact same explainer LLM and SAE that yielded the results of Section (ref) yielded the results in this section.

Conclusion

Literatures in empirical economics and econometrics have increasingly suggested the importance of open-ended discovery from high-dimensional or unstructured data (e.g., ludwig_machine-learning_2017,ludwig_machine_2024,chernozhukov_fisherschultz_2025). The framework proposed in this paper shows how newly validated statistical procedures for high-dimensional multiple hypothesis testing, when combined with the latest innovations in interpretability methods for machine learning models, can facilitate open-ended, interpretable discovery at scale with higher practicality and higher fidelity than previously possible.

The proposed framework also illustrates how treating (automatic) concept interpretation evaluation as a formal statistical inference problem can yield several scientific benefits in the process of making discoveries. A well-defined and intuitive quantitative method for evaluating the quality of concept descriptions is important in the context of the proposed framework, as it is these descriptions on which the scientific conclusions of discovery ultimately rest. By formulating autointerp evaluation as a statistical inference problem, the importance of evaluation with respect to the population of interest $P$ becomes salient; the need for assembling a held-out dataset becomes clear; the role for confidence intervals becomes legible by appeal to sampling uncertainty (c.f., miller_adding_2024); and we further unlock useful connections to the literature on post-selection inference, such that we may discipline settings evaluating a multiplicity of interpretations (c.f., andrews_inference_2024).

It is important to note that the proposed framework is compatible with many possible high-dimensional measurement technologies. The use of sparse dictionary learning models, and SAEs in particular, is viewed as a powerful default recommendation for creating high-dimensional concept embeddings: it is highly expressive and general as a measurement technology, and default adoption of such a technology further limits scope for motivated data mining. However, alternative choices of high-dimensional concept embedding may be more appealing depending on the use case of interest. An important future direction of work might also involve engineering a new “social science default” concept embedding that is tailored to a large space of solely social science relevant concepts.\footnote{We thank Jenny S. Wang for this suggestion.}

The proposed framework is most naturally viewed as one tool of many in the empirical researcher's toolkit for making discoveries from unstructured data. Using this framework alongside others that researchers may already be implementing is complementary, and would only serve to deepen insights into possible interpretations of inference on unstructured data.