Robust Offline Active Learning on Graphs

Consider a graph signal $\mathbf{f}\in\mathbb{R}^{n}$, where $\mathbf{f}(i)$ denotes the signal value at node $i$. For a set of nodes $S$, we define the subspace $\mathbf{L}_{\mathcal{S}}:=\{\mathbf{f}\in\mathbb{R}^{n}\mid\mathbf{f}(S^{c})=0\}$, where $\mathbf{f}(S)\in\mathbb{R}^{|\mathcal{S}|}$ represents the values of $\mathbf{f}$ on nodes in $\mathcal{S}$. In this paper, we consider both regression tasks, where $\mathbf{f}(i)$ is a continuous response, and classification tasks, where $\mathbf{f}(i)$ is a multi-class label.

Since $\mathbf{U}$ serves as a set of bases for $\mathbb{R}^{n}$, we can decompose $\mathbf{f}$ in the graph spectral domain as $\mathbf{f}=\sum_{j=1}^{n}\alpha_{\mathbf{f}}(\lambda_{j})U_{j}$, where $\alpha_{\mathbf{f}}(\lambda_{j})=\langle\mathbf{f},U_{j}\rangle$ is defined as the graph Fourier transform (GFT) coefficient corresponding to frequency $\lambda_{j}$. From a graph signal processing perspective, a smaller eigenvalue $\lambda_{k}$ indicates lower variation in the associated eigenvector $U_{k}$, reflecting smoother transitions between neighboring nodes. Therefore, the smoothness of $\mathbf{f}$ over the network can be characterized by the magnitude of $\alpha_{\mathbf{f}}(\lambda_{j})$ at each frequency $\lambda_{j}$. More formally, we measure the signal complexity of $\mathbf{f}$ using the bandwidth frequency $\omega_{\mathbf{f}}=\sup\{\lambda_{j}|\alpha_{\mathbf{f}}(\lambda_{j})>0\}$. Accordingly, we define the subspace of graph signals with a bandwidth frequency less than or equal to $\omega$ as $\mathbf{L}_{\omega}:=\{\mathbf{f}\in\mathbb{R}^{n}\mid\omega_{\mathbf{f}}\leq\omega\}$. It follows directly that $\forall\omega_{1}<\omega_{2},\mathbf{L}_{\omega_{1}}\subset\mathbf{L}_{\omega_{2}}$.

The key idea in graph-based semi-supervised learning is to reconstruct the graph signal $\mathbf{f}$ within a function space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ that depends on both the network structure and node-wise covariates, where the frequency parameter $\omega$ controls the size of the space to mitigate overftting. Assume that $Y_{i}$ is the observed noisy realization of the true signal $\mathbf{f}(i)$ at node $i$, active learning operates in a scenario where we have limited access to $Y_{i}$ on only a subset of nodes $\mathcal{S}$, with $|\mathcal{S}|<<n$. The objective is to estimate $\mathbf{f}$ within $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ using $\{Y_{i}\}_{i\in\mathcal{S}}$ by considering the empirical estimator of $\mathbf{f}$ as

where $l(\cdot)$ is a task-specific loss function. We denote $\mathbf{f}^{*}$ as the minimizer of (1) when responses on all nodes are available, i.e., $\mathbf{f}^{*}=\mathbf{f}_{\mathbf{V}}$. The goal of active semi-supervised learning is to design an appropriate function space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ and select an informative subset of nodes $\mathcal{S}$ for querying responses, under the query budget $|\mathcal{S}|\leq\mathcal{B}$, such that the estimation error is bounded as follows:

For a fixed $\mathcal{B}$, we wish to minimize the parameter $\rho>0$, which converges to $0$ as the query budget $\mathcal{B}$ approaches $n$.

In semi-supervised learning tasks on networks, both the network topology and node-wise covariates are crucial for inferring the graph signal. To effectively incorporate this information, we propose a function class for reconstructing the graph signal that lies at the intersection of the graph spectral domain and the space of node covariates. Motivated by the graph Fourier transform, we define the following function class:

Here, the choice of $\omega$ balances the information from node covariates and network structure. When $\omega=2$, $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ spans the full column space of covariates, i.e., $\text{Span}\{X_{1},\cdots,X_{p}\}$, allowing for a full utilization of the original covariate space to estimate the graph signal, but without incorporating any network information. On the other hand, when $\omega$ is close to zero—consider, for example, the extreme case where $|\{U_{j}\mid\lambda_{j}\leq\omega\}|=2$ and $p\gg 2$—then $\text{Proj}_{\mathbf{L}_{\omega}}X_{i}$ reduces to $\text{Span}\{U_{1},U_{2}\}$, resulting in a loss of critical information provided by the original $\mathbf{X}$.

By carefully choosing $\omega$, however, this function space can offer two key advantages for estimating the graph signal. From a signal recovery perspective, $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ imposes graph-based regularization over node covariates, enhancing generalizability when the dimension of covariates $p$ exceeds the query budget or even the network size—conditions commonly encountered in real applications. Additionally, covariate smoothing filters out signals in the covariates that are irrelevant to network-based prediction, thereby increasing robustness against potential noise in the covariates. From an active learning perspective, using $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ enables a bottom-up query strategy that begins with a small $\omega$ to capture the smoothest global trends in the graph signal. As the labeling budget increases, $\omega$ is adaptively increased to capture more complex graph signals within the larger space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$.

The graph signal $\mathbf{f}$ can be approximated by its projection onto the space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$. Specifically, let $\mathbf{U}_{d}=\{U_{1},U_{2},\cdots,U_{d}\}\in\mathbb{R}^{n\times d}$ stack the $d$ leading eigenvectors of $\mathcal{L}$, where $d=\operatorname*{argmax}_{1\leq j\leq n}(\lambda_{j}-\omega)\leq 0$. The graph signal estimation is then given by $\mathbf{U}_{d}\mathbf{U}^{T}_{d}\mathbf{X}{\beta}$, where $\beta\in\mathbb{R}^{d}$ is a trainable weight vector. However, the parameters $\beta$ may become unidentifiable when the covariate dimension $p$ exceeds $d$. To address this issue, we reparameterize the linear regression as follows:

where $\tilde{\beta}=\Sigma V_{2}^{T}{\beta}$ and $\tilde{\mathbf{X}}=\mathbf{U}_{d}V_{1}$. Here, $V_{1}\in\mathbb{R}^{d\times r}$, $V_{2}\in\mathbb{R}^{p\times r}$, and $\Sigma\in\mathbb{R}^{r\times r}$ denote the left and right singular vectors and the diagonal matrix of the $r$ singular values, respectively.

In the reparameterized form $(\ref{reg_2})$, the columns of $\tilde{\mathbf{X}}$ serve as bases for $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$, thus
$\text{dim}(\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A}))=\text{rank}(\tilde{\mathbf{X}})=r\leq\min\{d,p\}$. The transformed predictors $\tilde{\mathbf{X}}$ capture the components of the node covariates constrained within the low-frequency graph spectrum. A graph signal $\mathbf{f}\in\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ can be parameterized as a linear combination of the columns of $\tilde{\mathbf{X}}$, with the corresponding weights $\tilde{\beta}$ identified via

where $\tilde{\mathbf{X}}_{\mathbf{S}}\in\mathrm{R}^{|\mathcal{S}|\times r}$ is the submatrix of $\tilde{\bm{X}}$ containing rows indexed by the query set $\mathcal{S}$, and $\{\mathbf{f}(i)\}_{i\in\mathcal{S}}$ represents the true labels for nodes in $\mathcal{S}$. To achieve the identification of $\mathbf{f}$, it is necessary that $|\mathcal{S}|\geq r$; otherwise, there will be more parameters than equations in (3). More importantly, since $\text{rank}(\tilde{\mathbf{X}}_{\mathbf{S}})\leq\text{rank}(\tilde{\mathbf{X}})=r$, $\mathbf{f}$ is only identifiable if $\tilde{\mathbf{X}}_{\mathbf{S}}$ has full column rank. Notice that $r$ increases monotonically with $\omega$. If $\mathcal{S}$ is not carefully selected, the graph signal can only be identified in $\mathbf{H}_{\omega^{\prime}}(\mathbf{X},\mathbf{A})$ for some $\omega^{\prime}<\omega$, which is a subspace of $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$.

We first define the identification of $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ by the node query set $\mathcal{S}$ as follows:

Intuitively, the informativeness of a set $\mathcal{S}$ can be quantified by the frequency $\omega$ corresponding to the space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ that can be identified. To select informative nodes, we need to bridge the query set $\mathcal{S}$ in the spatial domain with $\omega$ in the spectral domain. To achieve this, we consider the counterpart of the function space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ in the spatial domain. Specifically, we introduce the projection space with respect to a subset of nodes $\mathcal{S}$ as follows: $\mathbf{H}_{\mathcal{S}}(\mathbf{X},\mathbf{A}):=\text{Span}\{{X}^{(\mathcal{S})}_{1},\cdots,{X}^{(\mathcal{S})}_{p}\}$, where

Here, ${X}_{p}(i)$ denotes the $p$-th covariate for node $i$ in $\mathbf{X}$. Theorem 3.1 establishes a connection between the two graph signal spaces $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ and $\mathbf{H}_{\mathcal{S}}(\mathbf{X},\mathbf{A})$, providing a metric for evaluating the informativeness of querying a subset of nodes on the graph.

We denote the quantity $\omega(\mathcal{S})$ in (4) as the bandwidth frequency with respect to the node set $\mathcal{S}$. This quantity can be explicitly calculated and measures the size of the space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ that can be recovered from the subset of nodes $\mathcal{S}$. The goal of the active learning strategy is to select $\mathcal{S}$ within a given budget to maximize the bandwidth frequency $\omega(\mathcal{S})$, thus enabling the identification of graph signals with the highest possible complexity.

To calculate the bandwidth frequency $\omega(\mathcal{S})$, consider any graph signal $\mathbf{g}$ and its components with non-zero frequency $\Lambda_{\mathbf{g}}:=\{\lambda_{i}\mid\alpha_{\mathbf{g}}(\lambda_{i})>0\}$. We use the fact that

where $\sum_{j:\lambda_{j}\in\Lambda_{\mathbf{g}}}c_{j}=1$ and $0\leq c_{j}\leq 1$. Combined with the Rayleigh quotient representation of eigenvalues, the bandwidth frequency $\omega_{\mathbf{g}}$ can be calculated as

As a result, we can approximate the bandwidth $\omega_{\mathbf{g}}$ using $\omega_{\mathbf{g}}(k)$ for a large $k$. Maximizing $\omega(\mathcal{S})$ over $\mathcal{S}$ then transforms into the following optimization problem:

where $\mathcal{B}$ represents the budget for querying labels. Due to the combinatorial complexity of directly solving optimization problem (5) by simultaneously selecting $\mathcal{S}$, we propose a greedy selection strategy as a continuous relaxation of (5).

The selection procedure starts with $\mathcal{S}=\emptyset$ and sequentially adds one node to $\mathcal{S}$ that maximizes the increase in $\omega(\mathcal{S})$ until the budget is reached. We introduce an $n$-dimensional vector $\mathbf{t}=(t_{1},t_{2},\cdots,t_{n})^{T}$ with $0\leq t_{i}\leq 1$, and define the corresponding diagonal matrix $D({\mathbf{t}})$ with diagonal entries given by $\mathbf{t}$. This allows us to encode the set of query nodes using $\mathbf{t}=\mathbf{1}_{\mathcal{S}}$, where $\mathbf{1}_{\mathcal{S}}(i)=1$ if $i\in\mathcal{S}$ and $\mathbf{1}_{\mathcal{S}}(i)=0$ if $i\in\mathcal{S}^{c}$. We then consider the space spanned by the columns of $D(\mathbf{t})\mathbf{X}$ as $\text{Span}\{D(\mathbf{t})\mathbf{X}\}$, and the following relation holds:

Intuitively, $\text{Span}\{D(\mathbf{t})\mathbf{X}\}$ acts as a differentiable relaxation of the subspace $\mathbf{H}_{\mathcal{S}}(\mathbf{X},\mathbf{A})$, enabling perturbation analysis of the bandwidth frequency when a new node is added to $\mathcal{S}$. The projection operator associated with $\text{Span}\{D(\mathbf{t})\mathbf{X}\}$ can be explicitly expressed as

To quantify the increase in $\hat{\omega}(\mathcal{S})$ when adding a new node to $\mathcal{S}$, we consider the following regularized optimization problem:

The penalty term on the right-hand side of (6) encourages the graph signal $\phi$ to remain in $\mathbf{H}_{\mathcal{S}^{c}}(\mathbf{X},\mathbf{A})$. As the parameter $\alpha$ approaches infinity and $\mathbf{t}=\mathbf{1}_{\mathcal{S}^{c}}$, the minimization $\lambda_{\alpha}(\mathbf{1}_{\mathcal{S}^{c}})$ in (6) converges to $\hat{\omega}(\mathcal{S})$ in (5). The information gain from labeling a node $i\in\mathcal{S}^{c}$ can then be quantified by the gradient of the bandwidth frequency as $t_{i}$ decreases from 1 to 0:

where $\phi$ is the minimizer of (6) at $\mathbf{t}=\mathbf{1}_{\mathcal{S}^{c}}$, which corresponds to the eigenvector associated with the smallest non-zero eigenvalue of the matrix $\mathbf{P}(\mathbf{1}_{\mathcal{S}^{c}})\mathcal{L}^{k}\mathbf{P}(\mathbf{1}_{\mathcal{S}^{c}})$. We then select the node $i=\operatorname*{argmax}_{j\in\mathcal{S}^{c}}\Delta_{j}$ and update the query set as $\mathcal{S}=\mathcal{S}\cup\{i\}$.

When calculating node-wise informativeness in (7), we can enhance computational efficiency by avoiding the inversion of $D(\mathbf{t})$. When $\mathbf{t}$ is in the neighborhood of $\mathbf{1}_{\mathcal{S}^{c}}$, we can approximate:

where $\mathbf{X}_{\mathcal{S}^{c}}=((X_{1})_{\mathcal{S}^{c}},\cdots,(X_{p})_{\mathcal{S}^{c}})$ and $\mathbf{Z}_{\mathcal{S}^{c}}=\mathbf{X}_{\mathcal{S}^{c}}(\mathbf{X}_{\mathcal{S}^{c}}^{T}\mathbf{X}_{\mathcal{S}^{c}})^{-1/2}$. Then, the node-wise informativeness can be explicitly expressed as:

We find that this approximation yields very similar empirical performance compared to the exact formulation in (7). Therefore, we adopt the formulation in (8) for the subsequent numerical experiments.

In real-world applications, we often have access only to a perturbed version of the true graph signals, denoted as $Y=\mathbf{f}+\mathbf{\xi}$, where $\mathbf{\xi}$ represents node labeling noise that is independent of the network data. When replacing the true label $\mathbf{f}(i)$ with $Y(i)$ in (3), this noise term introduces both finite-sample bias and variance in the estimation of the graph signal $\mathbf{f}$. As a result, we aim to query nodes that are sufficiently representative of the entire covariate space to bound the generalization error. To achieve this, we introduce principled randomness into the deterministic selection procedure described in Section 3.2 to ensure that $\mathcal{S}$ includes nodes that are both informative and representative. The modified graph signal estimation procedure is given by:

where $s_{i}$ is the weight associated with the probability of selecting node $i$ into $\mathcal{S}$.

Specifically, the generalization error of the estimator in (9) is determined by the smallest eigenvalue of $\tilde{\mathbf{X}}_{\mathcal{S}}^{T}\tilde{\mathbf{X}}_{\mathcal{S}}$, denoted as $\lambda_{\min}(\tilde{\mathbf{X}}_{\mathcal{S}}^{T}\tilde{\mathbf{X}}_{\mathcal{S}})$. Given that $\lambda_{\min}(\tilde{\mathbf{X}}^{T}\tilde{\mathbf{X}})=1$, our goal is to increase the representativeness of $\mathcal{S}$ such that $\lambda_{\min}(\tilde{\mathbf{X}}_{\mathcal{S}}^{T}\tilde{\mathbf{X}}_{\mathcal{S}})$ is lower-bounded by:

However, the informative selection method in Section 3.2 does not guarantee (10). To address this, we propose a sequential biased sampling approach that balances informative node selection with generalization error control.

The key idea to achieve a lower bound for $\lambda_{\min}(\tilde{\mathbf{X}}_{\mathcal{S}}^{T}\tilde{\mathbf{X}}_{\mathcal{S}})$ is to use spectral sparsification techniques for positive semi-definite matrices [4]. Let $v_{i}\in\mathbb{R}^{1\times r}$ denote the $i$-th row of the constrained basis $\tilde{\mathbf{X}}$. By definition of $\tilde{\mathbf{X}}$, it follows that $\mathbf{I}_{r\times r}=\sum_{i=1}^{n}v_{i}^{T}v_{i}$. Inspired by the randomized sampling approach in [21], we propose a biased sampling strategy to construct $\mathcal{S}$ with $|\mathcal{S}|\ll n$ and weights $\{s_{i}>0,i\in\mathcal{S}\}$ such that $\sum_{i\in\mathcal{S}}s_{i}v_{i}^{T}v_{i}\approx\mathbf{I}$. In other words, the weighted covariance matrix of the query set $\mathcal{S}$ satisfies $\lambda_{\min}(\tilde{\mathbf{X}}_{\mathcal{S}}^{T}W_{S}\tilde{\mathbf{X}}_{\mathcal{S}})\approx 1$, where $W_{S}$ is a diagonal matrix with $s_{i}$ on its diagonal.

We outline the detailed sampling procedure as follows. After the $(t-1)^{\text{th}}$ selection, let the set of query nodes be $\mathcal{S}_{t-1}$ with corresponding node-wise weights $\mathcal{W}_{t-1}=\{s_{j}>0\mid j\in\mathcal{S}_{t-1}\}$. The covariance matrix of $\mathcal{S}_{t-1}$ is given by $C_{t-1}\in\mathbb{R}^{r\times r}$, defined as $C_{t-1}=\tilde{\mathbf{X}}_{\mathcal{S}_{t-1}}^{T}\tilde{\mathbf{X}}_{\mathcal{S}_{t-1}}=\sum_{j\in\mathcal{S}_{t-1}}s_{j}v_{j}^{T}v_{j}$. To analyze the behavior of eigenvalues as the query set is updated, we follow [21] and introduce the potential function:

where $u_{t-1}$ and $l_{t-1}$ are constants such that $l_{t-1}<\lambda_{\min}(C_{t-1})\leq\lambda_{\max}(C_{t-1})<u_{t-1}$, and $\text{Tr}(\cdot)$ denotes the trace of a matrix. The potential function $\Phi_{t-1}$ quantifies the coherence among all eigenvalues of $C_{t-1}$. To construct the candidate set $B_{m}$, we sample node $i$ and update $C_{t}$, $u_{t}$, and $l_{t}$ such that all eigenvalues of $C_{t}$ remain within the interval $(l_{t},u_{t})$. To achieve this, we first calculate the node-wise probabilities $\{p_{i}\}_{i=1}^{n}$ as:

where $\sum_{i=1}^{n}p_{i}=1$. We then sample $m$ nodes into $B_{m}$ according to $\{p_{i}\}_{i=1}^{n}$. For each node $i\in B_{m}$, the corresponding weight is given by $s_{i}=\frac{\epsilon}{p_{i}\Phi_{t-1}},\;0<\epsilon<1$. After obtaining the candidate set $B_{m}$, we apply the informative node selection criterion $\Delta_{i}$ introduced in Section 3.2, i.e., selecting the node $i=\operatorname*{argmax}_{i\in B_{m}}\Delta_{i}$, and update the query set and weights as follows:

We then update the lower and upper eigenvalue bounds as follows:

The update rule ensures that $u_{t}-l_{t}$ increases at a slower rate than $u_{t}$, leading to the convergence of the gap between the largest and smallest eigenvalues of $\tilde{\mathbf{X}}_{\mathcal{S}}^{T}W_{S}\tilde{\mathbf{X}}_{\mathcal{S}}$, thereby controlling the condition number. Accordingly, the covariance matrix is updated with the selected node $i$ as:

With the covariance matrix update rule in (14), the average increment is $\mathbf{E}(C_{t})-C_{t-1}=\sum_{i=1}^{n}p_{i}s_{i}v_{i}^{T}v_{i}=\frac{\epsilon}{\Phi_{t-1}}\mathbf{I}$. Intuitively, the selected node allows all eigenvalues of $C_{t-1}$ to increase at the same rate on average. This ensures that $\lambda_{\min}(\tilde{\mathbf{X}}_{\mathcal{S}}^{T}\tilde{\mathbf{X}}_{\mathcal{S}})$ continues to approach $\lambda_{\min}(\tilde{\mathbf{X}}^{T}\tilde{\mathbf{X}})=1$ during the selection process, thus driving the smallest eigenvalue away from zero. Additionally, the selected node remains locally informative within the candidate set $B_{m}$. Compared with the entire set of nodes, selecting from a subset serves as a regularization on informativeness maximization, achieving a balance between informativeness and representativeness for node queries.

We summarize the biased sampling selection strategy in Algorithm 1. At a high level, each step in the biased sampling query strategy consists of two stages. First, we use randomized spectral sparsification to sample $m\ll n$ nodes and collect them into a candidate set $B_{m}$. Intuitively, the covariance matrix on the updated $\mathcal{S}$ maintains lower-bounded eigenvalues if a node from $B_{m}$ is added to $\mathcal{S}$. In the second stage, we select one node from $B_{m}$ based on the informativeness criterion in Section 3.2 to achieve a significant frequency increase in (7).

For initialization, the dimension of the network spectrum $d$, the size of the candidate set $m$, and the constant $0<\epsilon<1$ for spectral sparsification need to be specified. Based on the discussion at the end of Section 3.1, the dimension of the function space $\mathbf{H}_{\omega}(\mathbf{X},\mathbf{A})$ is at most $\mathcal{B}$, where $\mathcal{B}$ is the budget for label queries. Therefore, we can set $d=\min\{p,\mathcal{B}\}$. The parameters $m$ and $\epsilon$ jointly control the condition number $\frac{\lambda_{\max}(\tilde{\mathbf{X}}_{\mathcal{S}}^{T}W_{S}\tilde{\mathbf{X}}_{\mathcal{S}})}{\lambda_{\min}(\tilde{\mathbf{X}}_{\mathcal{S}}^{T}W_{S}\tilde{\mathbf{X}}_{\mathcal{S}})}$.

In practice, we can tune $m$ to ensure that the covariance matrix is well-conditioned. Specifically, we can run the biased sampling procedure multiple times with different values of $m$ and select the largest $m$ such that the condition number of the covariance matrix on the query set $\mathcal{S}$ is less than 10 [20]. This threshold is a commonly accepted rule of thumb for considering a covariance matrix to be well-conditioned [20]. Additionally, $\epsilon$ is typically fixed at a small value, following the protocol outlined in [21].

Based on the output from Algorithm 1, we solve the weighted least squares problem in (9):

and recover the graph signal on the entire network as $\hat{\mathbf{f}}=\tilde{\mathbf{X}}\hat{\beta}$.

Although our theoretical analysis is developed for node regression tasks, the proposed query strategy and graph signal recovery procedure are also applicable to classification tasks. Consider a $K$-class classification problem, where the response on each node $i$ is given by $\mathbf{f}(i)\in{1,2,\ldots,K}$. We introduce a dummy membership vector $(Y_{1}(i),\ldots,Y_{K}(i))$, where $Y_{c}(i)=1$ if $\mathbf{f}(i)=c$ and $Y_{c}(i)=0$ otherwise. For each class $c\in\{1,2,\ldots,K\}$, we first estimate $\hat{\beta}_{c}$ based on (15) with the training data $\{\tilde{\bm{X}}_{i\cdot},Y_{c}(i),s_{i}\}_{i\in\mathcal{S}}$, and then compute the score for class $c$ as $\hat{\mathbf{f}}_{c}=\tilde{X}\hat{\beta}_{c}$. The label of an unqueried node $j$ is assigned as $\mathbf{\hat{f}}(j)=\operatorname*{argmax}_{1\leq c\leq K}\{\hat{\mathbf{f}}_{1}(j),\hat{\mathbf{f}}_{2}(j),\cdots,\hat{\mathbf{f}}_{K}(j)\}$. Notice that the above score-based classifier is equivalent to the softmax classifier:

since the softmax function is monotonically increasing with respect to each score function $\{\mathbf{\hat{f}}_{c}\}_{c=1}^{K}$.

In the representative sampling stage, the computational complexity of calculating the sampling probability is $\mathcal{O}(n)$. We then sample $m$ nodes to formulate a candidate set $B_{m}$, where the complexity of sampling $m$ variables from a discrete probability distribution is $\mathcal{O}(m)$ [32]. Consequently, the complexity of the representative learning stage is $\mathcal{O}(n+m)$.

In the informative selection stage, we calculate the information gain $\Delta_{i}$ for each node. This involves obtaining the eigenvector corresponding to the smallest non-zero eigenvalue of the projected graph Laplacian matrix, with a complexity of $\mathcal{O}(n^{3})$ due to the singular value decomposition (SVD) operation. Subsequently, we compute $\Delta_{i}$ for each node in the candidate set $B_{m}$ based on their loadings on the eigenvector, which incurs an additional computational cost of $\mathcal{O}(mn)$. Therefore, the total complexity of our biased sampling method is $\mathcal{O}(n+m+nm+n^{3})$. Given a node label query budget $\mathcal{B}$, the overall computational cost becomes $\mathcal{O}(\mathcal{B}\left(n+m+nm+n^{3})\right)$.

When the dimension of node covariates $p\ll n$, we can replace the SVD operation with the Lanczos algorithm to accelerate the informative selection stage. The Lanczos algorithm is designed to efficiently obtain the $k$th largest or smallest eigenvalues and their corresponding eigenvectors using a generalized power iteration method, which has a time complexity of $\mathcal{O}(kn^{2})$ [17]. As a result, the complexity of the proposed biased sampling method reduces to $\mathcal{O}(pn^{2})$. This is comparable to GNN-based active learning methods, as GNNs and their variations generally have a complexity of $\mathcal{O}(pn^{2})$ per training update [5, 36].

We consider three different network topologies generated by widely studied statistical network models: the Watts–Strogatz model [34] for small-world properties, the Stochastic block model [16] for community structure, and the Barabási–Albert model [3] for scale-free properties.

Node responses are generated as $Y=\mathbf{f}+\xi$, where $\mathbf{f}$ is the true graph signal and $\xi\sim N(0,\sigma^{2}I_{n})$ is Gaussian noise. The true signal is constructed as $\mathbf{f}=\mathbf{U}_{d}\mathbb{\beta}$, where $\mathbb{\beta}$ is the linear coefficient and $\mathbf{U}_{d}$ denotes the leading $d$ eigenvectors of the normalized graph Laplacian of the synthetic network. Since our theoretical analysis assumes that the observed node covariates $\mathbf{X}$ contain noise, we generate $\mathbf{X}$ as a perturbed version of $\mathbf{U}_{d}$ by adding non-leading eigenvectors of the normalized graph Laplacian. The detailed simulation settings can be found in Appendix B.1.

We compare our algorithm with several offline active learning methods: 1. D-optimal [25] selects subset of nodes $\mathcal{S}$ to maximize determinant of observed covariate information matrix $\mathbf{X}_{\mathcal{S}}^{T}\mathbf{X}_{\mathcal{S}}$. 2. RIM [38] selects nodes to maximize the number of influenced nodes. 3. GPT [23] and SPA [13] split the graph into disjoint partitions and select informative nodes from each partition.