1 INTRODUCTION↩︎

The rapid advancement of network technologies has catalyzed the swift development and widespread deployment of intelligent tutoring systems (ITS), which aim to provide personalized and adaptive learning experiences [1]–[3]. Within the core framework of these systems, cognitive diagnosis (CD) serves as a fundamental component, dedicated to precisely modeling students’ latent knowledge states and skill proficiencies by analyzing their interactive learning behaviors and performance data.

Cognitive diagnosis methodologies can be broadly categorized into three main paradigms: psychometric theory-based methods [4]–[6], neural network-based methods [7]–[9], and graph modeling-based methods [10]–[12]. The first type, representing the foundational stage of research, primarily leverages psychological measurement techniques to assess learners’ latent abilities or traits. With the advent and flourishing of deep learning, researchers have increasingly turned to diverse neural network architectures to capture the complex, nonlinear student-exercise interaction patterns, thereby learning more effective representations. Furthermore, the critical importance of modeling both explicit and implicit high-order interactions among students, exercises, and knowledge concepts has gained significant recognition. This focus on relational complexity naturally motivates the introduction of graph structures into CD, leveraging their inherent strengths to model more comprehensive representations of the cognitive landscape.

Despite the remarkable progress achieved by existing CD methods, a critical limitation persists: they generally overlook the prevalent heterogeneous noise contamination in real educational environments. These approaches typically operate under the implicit assumption of relatively clean data, failing to explicitly address the adverse impacts of noisy data from various sources, including log miss, exercise confusion, \(Q\)-matrix confusion, and log flip, as illustrated in Figure 1. Furthermore, different model architectures exhibit varying levels of sensitivity to different noise sources, making it highly challenging to manually design a single, universally optimal structure capable of withstanding all types of heterogeneous noise contamination. As illustrated in Figure 2, we systematically modify the widely recognized NCD model by applying various replacements and ablations, resulting in multiple architectural variants. These variants exhibit significantly divergent performance change rates across diverse heterogeneous noise environments, directly validating that the noise tolerance of CD models is highly dependent on their underlying architecture.

To this end, we proposed OSCD, an evolutionary multi-objective One-Shot neural architecture search approach for Cognitive Diagnosis to explore efficient and robust model architectures. Specifically, OSCD operates through two sequential stages: training and searching. In the first stage, a search space encompassing diverse architectural combinations is constructed, and a weight-sharing supernet represented via the complete binary tree topology is trained. This process enables comprehensive exploration of potential architectures that transcend manual design priors. In the second stage, we explicitly formulate the search for the optimal architecture under heterogeneous noise scenarios as a MOP. To resolve this MOP, we develop an optimization framework that integrates a Pareto-optimal solution search strategy with cross-scenario performance evaluation. Finally, extensive experiments conducted on real-world educational datasets validate the effectiveness and robustness of the optimal architectures discovered by OSCD for CD tasks.

2 RELATED WORK↩︎

2.1 Cognitive Diagnosis↩︎

Cognitive diagnosis, as an advanced learner modeling approach, is theoretically grounded in the intersection of modern psychometrics and cognitive science. By analyzing students’ interactive learning behaviors, it can systematically deconstruct learners’ knowledge states and cognitive processes. Current research primarily focuses on manually developing sophisticated cognitive diagnostic models (CDMs) for more accurate learner modeling and assessment, which can be broadly classified into three categories: psychometric theory-based methods [4]–[6], neural network-based methods [7], [13], [14], and graph modeling-based methods [10], [15]–[17]. For instance, NCD[7], as a representative neural CDM, employs multidimensional parameters to construct fine-grained representations of students’ cognitive states and exercise attributes, while leveraging neural networks to model complex interactions. In addition, RCD[10] models student interactions and structural relationships by constructing a student-exercise-concept graph, employing graph convolutional networks (GCN[18]) to enhance the effectiveness of representation learning. More details of CD works are provided in the Appendix A.1.1.

With the increasing attention to learner modeling, recent years have also witnessed critical explorations from other perspectives, including interpretability[19], [20], robustness[21], and group learning[22], [23]. Despite the impressive results achieved by these cognitive diagnostic methods, their underlying CDM architectures often rely heavily on researchers’ domain expertise and experience. Such handcrafted models are limited in their ability to comprehensively explore all possible combinations of learning behaviors, thereby constraining the generalization and further investigation of cognitive diagnosis.

2.2 Neural Architecture Search↩︎

As a revolutionary automated machine learning technology[24], [25], Neural Architecture Search (NAS) has evolved into a comprehensive multi-domain and multi-modal framework since its pioneering works[26]–[28]. Currently, NAS techniques have achieved cross-domain penetration by offering automated optimization solutions for a wide range of mainstream deep neural networks (DNNs). This technology has been successfully applied in computer vision (CV)[29], natural language processing (NLP)[30], speech recognition systems[31], and graph neural networks (GNNs)[32]. Notably, the Transformer architecture, as a universal paradigm for cross-modal applications, has also benefited from the systematic exploration enabled by NAS methods[33]–[36]. More details of NAS works are provided in the Appendix A.1.2.

Meanwhile, several recent pioneering studies[37]–[39] have initiated the exploration of NAS techniques in the education domain. NASCD[40] adopts a NAS-based approach to search for CDMs, with its primary contribution lying in the design of a revised CD search space and the definition of an interpretability-oriented search objective, enabling the effective discovery of more interpretable and efficient models. Furthermore, ENAS-KT[41] achieves an optimal balance between input feature selection and architectural design by automatically selecting input features and formulating a novel Transformer-based search space. However, existing research has overly emphasized the high performance and interpretability of the searched CDMs, while insufficiently considering the robustness of these models in real-world scenarios. More critically, these methods have yet to incorporate computational complexity into their design considerations, leading to severe efficiency bottlenecks when scaling to complex interaction scenarios, which significantly hampers their developmental potential.

3 PRELIMINARY↩︎

3.1 Cognitive Diagnosis Task↩︎

In the field of intelligent education, let \(\mathcal{S}=\{s_1, s_2, \cdots, s_N\}\), \(\mathcal{E}=\{e_1, e_2, \cdots, e_M\}\), and \(\mathcal{C}=\{c_1, c_2, \cdots, c_K\}\) denote the sets of learners, exercises, and knowledge concepts, respectively, where \(N\), \(M\), and \(K\) represent the size of each set. The relationship between exercises and knowledge concepts is represented by a \(Q\)-matrix, denoted as \(\mathbf{Q}=\{q_{ij}\}^{M\times K}\), where \(q_{ij}=1\) if exercise \(e_i\) involves concept \(c_j\), and \(q_{ij}=0\) otherwise. The observed response data are provided in the form of the triplet set \(\mathcal{R}=\{(s_i,e_j,r_{ij})|s_i\in \mathcal{S},e_j\in \mathcal{E},r_{ij}\in \{0,1\}\}\), where \(r_{ij}=1\) denotes the learner \(s_i\) answers the exercise \(e_j\) correctly, and \(r_{ij}=0\) otherwise. Naturally, the response data \(\mathcal{R}\) is partitioned into training set \(\mathcal{R}_\mathrm{Tr}\), validation set \(\mathcal{R}_\mathrm{Val}\) and test set \(\mathcal{R}_\mathrm{Te}\). In this work, we construct multiple realistic heterogeneous noise scenarios. Specifically, we generate the perturbed validation sets \(\tilde{\mathcal{R}}_\mathrm{Val}=\{\mathcal{R}_\mathrm{Val}^\Theta, \mathcal{R}_\mathrm{Val}^\Phi, \mathcal{R}_\mathrm{Val}^\Psi, \mathcal{R}_\mathrm{Val}^\Omega \}\) by applying distinct noise perturbations \(\{ \Theta, \Phi, \Psi, \Omega \}\) to the original validation set \(\mathcal{R}_\mathrm{Val}\) derived from response data \(\mathcal{R}\). Among these, \(\Theta\) corresponds to partial response log miss, where \(\mathcal{R}_\mathrm{Val}^\Theta \subset \mathcal{R}_\mathrm{Val}\); \(\Phi\) reflects partial confusion for response exercises, where \(\mathcal{R}_\mathrm{Val}^\Phi=\{(s_i,e'_j,r_{ij})|e'_j\neq e_j,e'_j\in \mathcal{E}\}\); \(\Psi\) represents the confusion of the \(Q\)-matrix, where \(\mathbf{Q}=\{q'_{ij}\}^{M\times K}\) and \(q'_{ij}\neq q_{ij}\); \(\Omega\) denotes the flip of the response log, where \(\mathcal{R}_\mathrm{Val}^\Omega=\{(s_i,e_j,r'_{ij})|r'_{ij}\neq r_{ij},r'_{ij}\in \{0,1\}\}\). Typically, the noise interactions in real-world data remain implicit and unlabeled. Accordingly, the problem in this paper can be defined as follows:

Core cognitive diagnosis task. Given the learner-exercise response record \(\mathcal{R}\) and the predefined \(Q\)-matrix \(\mathbf{Q}\), the primary objective of cognitive diagnosis is to accurately and reliably estimate the learner’s proficiency on specific knowledge concepts.

Robust architecture search. Given the training set \(\mathcal{R}_\mathrm{Tr}\), the validation sets \(\mathcal{R}_\mathrm{Val} \cup \tilde{\mathcal{R}}_\mathrm{Val}\) and the test set \(\mathcal{R}_\mathrm{Te}\) of response records \(\mathcal{R}\) along with the \(Q\)-matrix \(\mathbf{Q}\), the primary objective of robust architecture search is to identify one or more model architectures that maintain stable and competitive performance across both original and perturbed validation sets. The searched architectures are then retrained from scratch on the training set and evaluated on the held-out test set \(\mathcal{R}_\mathrm{Te}\).

3.2 General Methodology of CD↩︎

For a general cognitive diagnostic task, the model accepts three main types of inputs and outputs a prediction score for student \(s_i\). In a general cognitive diagnosis framework, the model \(\mathcal{F}\) takes three primary types of inputs and outputs an estimated probability \(\hat{r}_{ij}\) for student \(s_i\) on exercise \(e_j\): \[\label{eq:1} \begin{gather} \hat{r}_{ij}=\mathcal{F}(\mathbf{h}_s,\mathbf{h}_e,\mathbf{h}_c), \\ \mathbf{h}_s = \mathbf{x}_i^s\times \mathrm{W}_s, \mathbf{h}_e = \mathbf{x}_j^e\times \mathrm{W}_e, \mathbf{h}_c = \mathbf{x}_j^e\times \mathrm{Q},\\ \end{gather}\tag{1}\] where \(\mathbf{h}_s \in \mathbb{R}^{1\times d}\) is the learner-related feature embedding, \(\mathbf{h}_e \in \mathbb{R}^{1\times d}\) is the exercise-related feature embedding, \(\mathbf{h}_c \in \mathbb{R}^{1\times K}\) is the concept-related feature embedding, \(d\) denotes the embedding dimension, \(\mathbf{x}_i^s \in \{0,1\}^{1\times N}\) denotes the one-hot vector of student \(s_i\), \(\mathbf{x}_i^e \in \{0,1\}^{1\times M}\) denotes the one-hot vector of exercise \(e_j\), \(\mathrm{W}_s\) and \(\mathrm{W}_e\) are trainable matrices.

4 METHODOLOGY↩︎

4.1 Overview of The Proposed OSCD↩︎

As illustrated in Figure 3 (see Algorithm 1 for stepwise summary in Appendix A.2), our proposed OSCD model comprises two core stages: training and searching. During the training stage, we first design a search space for CD tasks consisting of two components: the first integrates three types of input nodes and \(17\) internal computational nodes, while the second dynamically selects the final prediction layer based on the output dimension from the preceding part. We then construct and train a supernet based on this search space. To fully explore the combinatorial potential of models, we refine the search space into a complete binary tree topology representation. This structured representation significantly facilitates the development and evaluation of all possible model combinations.

Upon completion of supernet training, the search stage commences by initializing a population of \(Pop\) individuals. Some individuals are initialized using existing high-performance manually-designed architectures to leverage their empirical value and accelerate search convergence, while the remaining individuals are randomly generated from the search space. Subsequently, parent individuals are selected through standard binary tournament selection to form a mating pool, followed by the application of single-point crossover and bit-wise mutation to generate offspring populations. Using the trained supernet weights, we efficiently evaluate the performance of all newly generated offspring architectures. Then, employing the environmental selection mechanism of NSGA-II [42], we filter individuals with superior fitness from the combined parent and offspring populations to form a new generation. This search cycle iterates until reaching the maximal generation count, ultimately outputting all non-dominated individuals in the final population as the Pareto-optimal solution set.

4.2 Training Stage↩︎

4.2.1 Search space↩︎

As mentioned earlier, the search spaces employed by existing NAS methods are task-specific and cannot be directly applied to the CD task[25], [27].

Therefore, we systematically investigate and summarize the architectural characteristics of existing CDMs. They can be uniformly abstracted into a general paradigm comprising multiple input nodes, multiple internal nodes, and one output node, as illustrated in Figure 4 (a).

To enable the general model to represent as many CDMs as possible, we systematically integrated \(17\) categories of representative feature transformation operations, which are structurally embedded in the model’s internal computational nodes to support general representations of diverse CDMs. Table 1 presents the details of these \(17\) operators and their compliance with the Lipschitz condition, which are primarily categorized into three types: a placeholder operation(applicable to both input and operation), unary operations (one input), and binary operations (two inputs). It should be specially noted that we have categorized five operators: \(Sum\), \(Mean\), \(FFN\), \(FFN\_D\), and \(Concat\), collectively termed as dimension-sensitive operators. The key characteristic of these operators is that they require \(d\)-dimension vectors as input to execute their specific computational logic. Specifically, to address the critical challenge of potential partial path missing in tree-like topological structures within the search space, we innovatively introduce a placeholder operation (\(Zero\)) mechanism. This mechanism maintains the continuity of feature transmission in incomplete binary tree structures by constructing virtual paths, ultimately enabling complete cognitive state derivation from leaf nodes to the root node. Furthermore, we conduct rigorous verification of Lipschitz continuity for these \(17\) operator categories, precisely determining the upper bounds of their Lipschitz constants. This mathematical property directly characterizes the output stability of nodes under heterogeneous noise perturbation. Additional details regarding the operators and proofs of the Lipschitz condition are provided in the Appendix A.3 due to space constraints.

As shown in Table 1, all operations output either a scalar or a \(d\)-dimension vector, indicating that the general model must ultimately produce either a scalar \(y\) or a \(d\)-dimension vector \(\mathbf{y}\). To accomplish the prediction task, we need to specify corresponding computational processes for both scalar and \(d\)-dimension vector outputs, as detailed below: \[\hat{r}_{ij}=\left\{ \begin{align} & y,~if~y\in \mathbb{R}^1 \\ \mathrm{FC}_3( & \mathrm{FC}_2(\mathrm{FC}_1(\mathbf{y})),~if~y\in \mathbb{R}^{1\times d} \end{align} \right.,\] where \(\mathrm{FC}_1(\cdot)\), \(\mathrm{FC}_2(\cdot)\), and \(\mathrm{FC}_3(\cdot)\) are three Fully Connected(FC) layers with output dimensions \(512\), \(256\), \(1\), respectively. It can also maintain certain monotonicity properties, which has gained recognition from numerous researchers and has been widely adopted in most CD studies[7], [10].

Therefore, the variable components in the general model constitute the search space for CD, primarily consisting of two parts. The first part mainly includes input nodes and internal computational nodes that form the primary computational path. The second part employs an identity operation and a three-layer FC network to accomplish the prediction task. As illustrated in Figure 4 (a), the first part contains up to three types of input nodes (represented by rounded rectangles) and \(17\) internal computational nodes (represented by ovals, with red indicating binary operations and yellow indicating unary operations). The second part primarily consists of an identity operation (represented by a gray rectangle) and a three-layer FC network (represented by yellow rectangles) to obtain the final predicted value \(\hat{r}_{ij}\).

4.2.2 Supernet↩︎

However, while the search space incorporating diverse input nodes and internal computational nodes offers expressive flexibility, it suffers from a critical flaw: the combinatorial explosion creates an excessively large search space that imposes prohibitively high computational complexity on traditional NAS methods. Specifically, the independent training and evaluation processes for numerous substructures are computationally intensive, and the discrete search mechanism struggles to efficiently converge to high-performance regions.

To overcome computational efficiency bottlenecks and accelerate optimization convergence, we innovatively introduce a supernet architecture. Constructed upon the aforementioned search space, this architecture employs weight-sharing mechanisms and layered sandwich training rules to achieve synchronous representation and evaluation of massive subnets, substantially reducing the time complexity of traditional NAS methods and establishing an efficient computational foundation for the subsequent MOP. Based on the above, we employ a tree-based topological structure to encode CDMs, as shown in Figure 4 (b). Specifically, we construct a full binary tree topology with height \(h\), where all input nodes serve as leaf nodes, intermediate nodes function as non-leaf nodes, and the output node corresponds to the root node. This results in a full binary tree topology consisting of \(2^{h-1}\) leaf nodes and \(2^{h-1}-1\) non-leaf nodes, representing \(2^{h-1}\) input nodes (\(Zero\), \(\mathbf{h}_s\),\(\mathbf{h}_e\), and \(\mathbf{h}_c\)) and \(2^{h-1}-1\) operation nodes (\(17\) operators in Table 1), respectively. Naturally, for each architecture, its \(i\)-th operation node can be represented as \(node_i=\{input_i^L, input_i^R, operator_i\}\), where \(input_i^L\), \(input_i^R\), and \(operator_i\) respectively denote the computational result of the left subtree, the right subtree, and the currently selected operator. The prediction result \(\hat{r}_{ij}\) is ultimately generated through recursive computation of the feature propagation process from leaf nodes to the root node.

However, reality proves less ideal as this full binary tree representation paradigm encounters significant challenges when generalizing to complex CDMs: (1) Dimension Incompatibility: This primarily manifests as matching errors between input nodes and operators, which can be categorized into four typical scenarios: (a) Unary operator mismatch: as illustrated in Figure 5 (a), when a unary operation (\(Sigmoid\)) is incorrectly connected to two inputs \(( \mathbf{h}_s, \mathbf{h}_e)\), we enforce computation using the left input: \(sigmoid(\mathbf{h_s})\). (b) Binary operator deficiency: as shown in Figure 5 (b), when a binary operation (\(Concat\)) receives only a valid left input \(\mathbf{h}_s\) while the right input is an invalid placeholder \(Zero\), we pad the right dimension with an all-ones tensor: \(Concat(\mathbf{h}_s,\mathbf{1}^{1\times d})\). (c) Dimension-sensitive operator left input error: as depicted in Figure 5 (c), when a dimension-sensitive operator (\(FFN\)) receives a scalar \(Sum(\mathbf{h}_s)\) as left input but requires a \(d\)-dimensional vector, we directly substitute it with a dimensionless operator to avoid dimensional conflicts: \(Sigmoid(Sum(\mathbf{h}_s))\). (d) Dimension-sensitive operator right input error: as can be seen from Figure 5 (d), when the left input \(Sigmoid(\mathbf{h}_s)\in \mathbb{R}^{1\times d}\) but the right input is a scalar \(Sum(\mathbf{h}_e)\in \mathbb{R}^1\) (in \(Concat\) operations), we expand the right input dimension using an all-ones tensor: \(Concat(Sigmoid(\mathbf{h}_s), [Sum(\mathbf{h}_e)|| \mathbf{1}^{1\times (d-1)}])\). Detailed explanations are provided in Appendix A.4.

(2) Path Discontinuity: A significant portion of computational paths originating from leaf nodes fail to reach the root node, as shown in Figure 5 (e). This occurs because the computational graph contains numerous \(Zero\) placeholders acting as operators. To address this, we preserve the active computational flows along the paths. \(Zero\) placeholders are tagged accordingly, and these tags are subsequently used to validate path viability and determine the final output.

(3) Structural Collapse: As illustrated in Figure 5 (f), when computational paths cannot extend to the root node, or when both the root node and all its child nodes are \(Zero\) placeholders, the model loses its effective predictive capability. To address structural collapse that disables prediction functionality, we introduce a dynamic fallback mechanism. When it detects either that no computation path originating from leaf nodes can reach the root, or the root node itself is a \(Zero\) placeholder, it automatically activates a contingency procedure. This procedure traces back to the highest-level valid computational node proximal to the root and adopts its output as the final prediction. This mechanism elevates system fault tolerance by ensuring graceful degradation under structural collapse, thereby preserving continuous prediction capability without compromising the core computational paradigm.

Building upon the solutions to the aforementioned issues and inspired by the sandwich rule [29], we propose employing three types of subnets in each training batch of the supernet: the minimal subnet (with the smallest parameters representing the lower bound), the maximal subnet (with the largest parameters representing the upper bound), and multiple randomly sampled subnets (typically \(3\) samples). Consequently, this approach not only emphasizes training boundaries but also introduces randomness during supernet training, achieving a balance between extreme search and fair exploration of the search space, thereby enabling the trained supernet weights to better facilitate subsequent subnet performance evaluation. The Appendix A.5 presents the specific details.

4.3 Searching Stage↩︎

4.3.1 Individual Representation and Population Initialization↩︎

Following prior works, we first define individual representation. To precisely encode candidate architectures from the search space (Figure 4) into evolvable individuals, we utilize a fixed-length integer vector \(\mathbf{Individual}_i=[ind_1,\dots,ind_{2^{h-1}},ind_{2^{h-1}+1},\dots,ind_{2^h-1}]\) where \(h\) is the height of the full binary tree topology. Each element in \(\mathbf{Individual}_i\) corresponds to a configurable decision point within the search space, such as input selection or operator choice. The integer value at each position explicitly specifies the concrete configuration option. Most crucially, we design a tailored population initialization mechanism. When constructing the initial population of size \(Pop\), we concurrently deploy: uniform random sampling for comprehensive search space coverage, and injection of elite genes from validated cognitive diagnosis models such as IRT, MIRT, and NCD. This dual-strategy framework leverages prior knowledge to accelerate convergence while maintaining topological diversity, ultimately enhancing Pareto frontier approximation efficiency without compromising global search capacity.

4.3.2 Objective Function↩︎

To simultaneously optimize candidate architectures for accuracy on original validation data and robustness in complex heterogeneous noise scenarios while ensuring the distribution characteristics of model predictions under noise interference closely match those on original data, we define the following three objective functions for joint optimization to enhance behavioral consistency and stability at the output level: \[\max_{Archi}F(Archi)=\left\{ \begin{align} f_1(Archi) & =AUC(Archi; \mathcal{R}_\mathrm{Val}) \\ f_2(Archi) & =AUC(Archi; \tilde{\mathcal{R}}_\mathrm{Val}) \\ f_3(Archi) & =D_{\mathrm{KL}}(\mathcal{P} || \tilde{\mathcal{P}};Archi) \end{align} \right.,\] where \(Archi\) represents the candidate architecture to be optimized, \(AUC(Archi;D)\) represents the AUC (Area Under an ROC Curve) of \(Archi\) on the dataset \(D\), \(\mathcal{P}\) and \(\tilde{\mathcal{P}}\) represents the predicted probability distributions of \(archi\) on \(\mathcal{R}_\mathrm{Val}\) and \(\tilde{\mathcal{R}}_\mathrm{Val}\) respectively, and \(D_{\mathrm{KL}}\) is the Kullback-Leibler divergence between these two probability distributions. This objective function drives the search process to simultaneously optimize three critical dimensions:

(1) Fundamental Performance Guarantee: \(f_1\) uses AUC as the performance baseline to ensure the model possesses essential discriminative characteristics in clean data domains. This dimension directly reflects the expected prediction accuracy of the model in ideal environments.

(2) Noise Robustness Characterization: \(f_2\) systematically quantifies operational fluctuations in non-ideal environments by constructing heterogeneous noise scenarios \(\tilde{\mathcal{R}}_\mathrm{Val}=\{\mathcal{R}_\mathrm{Val}^\Theta, \mathcal{R}_\mathrm{Val}^\Phi, \mathcal{R}_\mathrm{Val}^\Psi, \\ \mathcal{R}_\mathrm{Val}^\Omega \}\) (including multi-source corruptions like Log Miss, Exercise Confusion, Q-matrix Confusion, and Log Flip). This reflects the model’s performance instability when facing heterogeneous perturbations.

(3) Distribution Consistency Constraint: \(f_3\) employs KL divergence minimization as the distribution alignment criterion, enforcing functional convergence between clean data prediction perturbation \(\mathcal{P}\) and noisy environment perturbation \(\tilde{\mathcal{P}}\). This mechanism explicitly penalizes prediction behavior deviations, fundamentally ensuring the topological invariance of model behavior under perturbations. Specifically, by minimizing the KL divergence between distributions, we effectively encourage the model to maintain consistency in its output distribution when exposed to perturbations or noise. This mechanism prevents the model from overfitting to random fluctuations in individual samples, thereby mitigating the risk of reduced robustness. The proof is provided in Appendix A.6.

4.3.3 Genetic Operations↩︎

For fixed-length integer vectors (length \(2^h-1\)) representing full binary trees of height h, we employ the standard single-point crossover and the bit-wise mutation [43] operations. The single-point crossover randomly selects a position in the vector and exchanges the gene segments before and after this position between two parent individuals to produce offspring. Then, the bit-wise mutation randomly selects a single gene locus in the vector and modifies its value to another valid value within the defined domain. These operations, performed within the closed integer vector space \(\mathbb{Z}^{2^h-1}\), maintain an isomorphic mapping to the tree structure through the encoding scheme, strictly preserving the topological integrity of the full binary tree.

5 EXPERIMENT↩︎

In this section, we conduct comprehensive experiments on two publicly available real-world educational datasets to validate the effectiveness of our proposed approach. Specifically, we will answer the following research questions (RQs) to unfold the experiments:

• \(\boldsymbol{RQ1}\): How does the proposed OSCD model demonstrate effectiveness and superiority in addressing the CD task?

• \(\boldsymbol{RQ2}\): How does the robustness of OSCD compare to other baseline methods under different types of noise?

• \(\boldsymbol{RQ3}\): Can OSCD’s supernet reliably predict the relative performance of random architectures?

• \(\boldsymbol{RQ4}\): Is OSCD efficient in terms of architecture search time?

• \(\boldsymbol{RQ5}\): What insights can be gained from visualizing OSCD’s searched architectures and their objective function values?

5.1 Experimental Setting↩︎

5.1.1 Datasets.↩︎

To assess the effectiveness of the proposed model in learner modeling, we conducted experiments on two publicly available real-world educational datasets: ASSIST09[44], SLP-Math[45]. More details of the dataset descriptions are available in Appendix A.8.

5.1.2 Baseline Methods.↩︎

We conducted extensive experiments and compared several strong and commonly used baselines, including DINA[4], IRT[5], MIRT[6], MF[46], NCD[7], RCD[10], KSCD[13], KaNCD[47], ReliCD[48], ORCDF[16], and EMO-NAS-CD[40]. The introduction and implementation details can be found in Appendix A.9. We evaluated the models utilizing three common metrics, including ACC, AUC, and RMSE.

5.1.3 Implementation Details.↩︎

Following the approach described in [7], [49], for each dataset, we filtered out learners with fewer than \(15\) historical response records to ensure that every learner has sufficient exercise records for diagnosis. All experimental results are reported as averages computed over ten independent runs. We implemented all methods with PyTorch by Python, and the experiments were conducted on an NVIDIA Geforce RTX4090 GPU. For the baselines, we utilized the optimal parameter configurations suggested in the original work. In training stage, the embedding size \(d\), batch size, and learning rate were set to \(128\), \(128\), \(0.001\), respectively. The height \(h\) of tree was set to \(6\). In searching stage, the population size \(Pop\) was set to \(100\), the maximal number of generation was set to \(100\). Our codes are available at https://github.com/AhuTwelve/OSCD.

5.2 Performance Comparison↩︎

Table ¿tbl:tab:main95experiment? compares the ACC, AUC, and RMSE of the baselines, the top-performing architectures discovered by EMO-NAS-CD, and those discovered by our proposed OSCD on two educational datasets: ASSIST09 and SLP-Math. The best results for each model are highlighted in bold, while the second-best results are underlined. Table ¿tbl:tab:main95experiment? (a) and Table ¿tbl:tab:main95experiment? (b) respectively present the performance comparisons of various methods under different datasets and their corresponding split ratios, aiming to validate the robustness and effectiveness of our method across diverse experimental settings. According to the results, there are several observations: (1) The architectures discovered by our proposed OSCD method achieved the best or second-best performance on both datasets, significantly outperforming most existing baseline models overall, which fully demonstrates the effectiveness of the proposed method. Notably, on the SLP-Math dataset, the discovered architectures SA1 and SA5 achieved the best results across almost all evaluation metrics, further validating OSCD’s capability in architecture optimization when dealing with real-world educational data. (2) Under various split settings, the discovered architectures maintained stable performance across different metrics, showing minimal fluctuations and demonstrating strong generalization and resistance to overfitting. Compared with certain baselines such as IRT, MF, and KaNCD, which suffer from performance degradation under limited training data, the robustness of the models discovered by OSCD is particularly noteworthy.

5.3 Robustness Study↩︎

To verify the robustness of the model architectures obtained by OSCD, we simulate different heterogeneous noise environments on the original test set and inject noise at varying ratios. By comparing the AUC change ratios under different noise ratios, we can effectively assess the model’s stability under noise interference. A positive change ratio indicates that the model’s AUC improves under noise interference, while a negative change ratio indicates a decline in AUC. The absolute value of the change ratio reflects the magnitude of performance fluctuations under noise; the smaller the absolute value, the less sensitive the model is to noise, and hence, the stronger its robustness. As illustrated in Table 2, there are several noteworthy observations: (1) Across all noise types and ratio settings, the architectures SA1 to SA5 discovered by OSCD consistently exhibit the smallest average AUC change ratios, indicating the lowest susceptibility to noise and the smallest performance fluctuations. Compared with other baseline methods, OSCD achieves significantly lower fluctuation magnitudes in most cases, demonstrating stronger robustness and stability. (2) Regarding the impact of heterogeneous noise types, log miss and \(Q\)-matrix confusion have relatively minor interference with model performance, whereas the other two types of noise cause severe damage to the models, particularly log flip, which can degrade the predictive performance of most models to a level close to random guessing.

5.4 Ranking Correlation Evaluation↩︎

To evaluate the fidelity of the supernet’s performance prediction, we randomly generated 50 architectures and compared the rankings obtained from direct supernet evaluation with those from independently trained models, calculating Kendall’s tau and Spearman’s rank correlation coefficients to quantify the consistency between the two evaluation results. As depicted in the Figure 6, it can be observed that the performance rankings predicted by the supernet are highly consistent with those obtained from independent training across the ACC, AUC, and RMSE evaluation metrics, with most data points lying along the diagonal, indicating strong ranking consistency. Moreover, the Kendall’s tau and Spearman’s rank correlation coefficients of AUC are close to \(1\), further confirming the high fidelity of the supernet in performance prediction. This demonstrates that OSCD’s supernet can reliably reflect the true relative performance of random architectures, providing a solid foundation for efficient and accurate architecture search.

5.5 Time Efficiency Analysis↩︎

Table 3 presents the runtime comparison on ASSIST09 and SLP-Math. While traditional CDMs run quickly, their limited expressiveness constrains performance. EMO-NAS-CD achieves strong results but demands \(15.2\) GPU days on ASSIST09 and \(30.6\) hours on SLP-Math, making large-scale exploration impractical. OSCD delivers a decisive efficiency breakthrough. It reduces runtime to \(26.1\) hours on ASSIST09 (\(14×\) faster) and \(51.7\) minutes on SLP-Math (\(35×\) faster), enabling rapid architecture exploration within feasible time budgets. This efficiency transforms NAS in cognitive diagnosis from a costly niche method into a practical, scalable tool, unlocking more iterative experimentation and broader research opportunities. By lowering the computational barrier, OSCD paves the way for more agile, exploratory, and impactful research in robust learner modeling and intelligent tutoring systems.

5.6 Visualization of Discovered Architectures↩︎

To further investigate the reasons behind the robustness of the searched architectures, we conducted additional visualization experiments. As shown in the Figure 7 and 8, the visualization results reveal that the non-dominated individuals generally select operations that satisfy the Lipschitz condition. Such operations can effectively constrain the sensitivity of the network outputs to input perturbations, which is a key factor in enhancing robustness. This phenomenon offers a new perspective for understanding the formation mechanism of robust architectures and may provide further inspiration for the design of future robust neural networks.

6 CONCLUSION↩︎

In this work, we proposed OSCD, an evolutionary multi-objective One-Shot neural architecture search framework tailored for Cognitive Diagnosis, aiming to identify both efficient and robust model architectures. Specifically, OSCD proceeds in two stages: a training stage and a searching stage. During training, a comprehensive search space covering diverse architectural configurations is established, and a weight-sharing supernet based on a complete binary tree topology is optimized, enabling the exploration of architectures beyond hand-crafted design constraints. In the subsequent searching stage, the problem of identifying the optimal architecture under heterogeneous noise conditions is formulated as a multi-objective optimization task. This is addressed through a framework that combines Pareto-optimal solution discovery with cross-scenario performance evaluation. Extensive experiments on real-world educational datasets demonstrate that the architectures obtained through OSCD consistently deliver effective and robust performance in CD tasks, underscoring the practical value of the proposed approach. We hope this work will provide new perspectives in the field of CD.

7 Appendix↩︎

7.1 RELATED WORK↩︎

7.1.1 Cognitive Diagnosis↩︎

Cognitive diagnosis, as an advanced learner modeling approach, is theoretically grounded in the intersection of modern psychometrics and cognitive science. By analyzing students’ interactive learning behaviors, it can systematically deconstruct learners’ knowledge states and cognitive processes. Current research primarily focuses on manually developing sophisticated cognitive diagnostic models (CDMs) for more accurate learner modeling and assessment, which can be broadly classified into three categories: psychometric theory-based methods, neural network-based methods, and graph modeling-based methods. In the early stages of cognitive diagnosis development, researchers mainly employed psychometric techniques to assess learners’ abilities, utilizing unidimensional or multidimensional latent vectors to measure skill proficiency, as exemplified by DINA[4], IRT[5], and MIRT[6]. The advancement of deep learning has sparked significant interest, leading to the emergence of various neural network-based approaches [7], [13], [14]. For instance, NCD[7], as a representative neural CDM, employs multidimensional parameters to construct fine-grained representations of students’ cognitive states and exercise attributes, while leveraging neural networks to model complex interactions. KSCD[13] feeds student, exercise, and concept vectors into an IRT-like neural network as a predictive diagnostic function to uncover implicit relationships among concepts. Similarly, domain-level zero-shot cognitive diagnosis[14] and diagnosis using incremental response logs[50] have been developed to address practical challenges in online learning environments.

Owing to the unique characteristics of graph structures and the advancements in graph learning, graph modeling-based methods[10], [15]–[17] have gained widespread acceptance. These methods primarily leverage graph topological information to reveal higher-order relationships among students, exercises, and concepts, as well as to mine valuable latent implicit associations. For example, RCD[10] models student interactions and structural relationships by constructing a student-exercise-concept graph, employing graph convolutional networks (GCN[18]) to enhance the effectiveness of representation learning. With the increasing attention to learner modeling, recent years have also witnessed critical explorations from other perspectives, including interpretability[19], [20], robustness[21], and group learning[22], [23].

7.1.2 Neural Architecture Search↩︎

As a revolutionary automated machine learning technology[24], [25], Neural Architecture Search (NAS) has evolved into a comprehensive multi-domain and multi-modal framework since its pioneering works[26]–[28]. Currently, NAS techniques have achieved cross-domain penetration by offering automated optimization solutions for a wide range of mainstream deep neural networks (DNNs). This technology has not only been successfully applied to the construction of core models like CNNs in computer vision (CV)[29], but has also extended to recurrent neural networks (RNNs) in natural language processing (NLP)[30], speech recognition systems[31], and graph neural networks (GNNs) for tasks such as node classification and link prediction[32]. Notably, the Transformer architecture, as a universal paradigm for cross-modal applications, has also benefited from the systematic exploration enabled by NAS methods[33]–[36].

7.2 Pseudocode of OSCD↩︎

Next, we present the core algorithm pseudocode of OSCD to facilitate a thorough understanding of its internal working mechanism.

7.3 Proof of Lipschitz Condition↩︎

We provide more detailed supplementary explanations for the operators in the table. First, it is worth noting that \(FFN\), \(FFN\_D\) and \(Concat\) are all FC-layer-based operations, indicating they incorporate learnable parameters to enhance the representational capacity of the model architecture. Moreover, while other operators contain no learnable parameters, the embedding layer in Eq.@eq:eq:1 ensures the general model avoids being limited to simple symbolic computations. Next, we introduce the definition and proof of the Lipschitz condition.

Definition 4.1(Lipschitz Condition). Let \(f: D\rightarrow \mathbb{R}\) be a function defined on a set \(D\subseteq \mathbb{R}\). The function \(f\) is said to satisfy that Lipschitz condition (or be Lipschitz continuous) on \(D\) if there exists a constant \(L \geq 0\) such that \[|f(x)-f(y)| \leq L|x-y|, \forall x,y\in D.\] Any constant \(L\) satisfying this inequality is called a Lipschitz constant for \(f\) on \(D\). The Lipschitz constant of \(f\), denoted \(L_f\), is the infimum (greatest lower bound) of all Lipschitz constants for \(f\) on \(D\); that is, \[L_f:=inf\{L \geq 0 ||f(x)-f(y) \leq L|x-y|,\forall x,y\in D\}.\]

Proof:

• Iden(): Given that \(x\), \(y\) and \(Iden(x)=x\), we then have \[|Iden(x)-Iden(y)|=|x-y|\leq 1\cdot |x-y|.\] Therefore, it can be proved that \(Iden\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=1\).

• Neg(): Given that \(x\), \(y\) and \(Neg(x)=-x\), we then have \[|Neg(x)-Neg(y)|=|y-x|=|x-y|\leq 1\cdot |x-y|.\] Therefore, it can be proved that \(Neg\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=1\).

• Abs(): Given that \(x\), \(y\) and \(Abs(x)=|x|\), we then have \[|Abs(x)-Abs(y)|=||x|-|y||\leq |x-y|= 1\cdot |x-y|.\] Therefore, it can be proved that \(Abs\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=1\).

• Inv(): Given that \(x\), and \(Inv(x)=1/(x+1e^{-6})\), we can obtain \[Inv'(x)=-\frac{1}{(x+1e^{-6})^2},\] when \(x\rightarrow -1e^{-6}\), \(|Inv'(x)|\rightarrow \infty\). Based on the definition of Lipschitz condition, if the derivative of a function is unbounded, then the function does not satisfy the Lipschitz condition. Therefore, it can be proved that \(Inv\) operator does not satisfy the Lipschitz condition.

• Square(): Given that \(x\), and \(Square(x)=x^2\), we can obtain \[Square'(x)=2x,\] when \(x\rightarrow \infty\), \(|Square'(x)|\rightarrow \infty\). According to the definition of Lipschitz condition, if the derivative of a function is unbounded, then the function does not satisfy the Lipschitz condition. Therefore, it can be proved that \(Square\) operator does not satisfy the Lipschitz condition.

• Sqrt(): Given that \(x\), and \(Sqrt(x)=sign(x)\cdot \sqrt{|x|+1e^{-6}}\), we can obtain \[Sqrt'(x)=\left\{ \begin{align} \frac{1}{2\sqrt{x+1e^{-6}}}, x\geq 0 \\ \frac{1}{2\sqrt{1e^{-6}-x}}, x<0 \end{align} \right.,\] when \(x\rightarrow 0\), \(Sqrt'(x)\rightarrow 1/{2\sqrt{1e^{-6}}}\), where \(1e-6\) is the special value we set. Theoretically, \(Sqrt(x)\) is not continuous at \(x=0\) in terms of its derivative. Therefore, we believe that \(Sqrt\) operator does not satisfy the Lipschitz condition.

• Tanh(): Given that \(x\), \(y\) and \(Tanh(x)=tanh(x)=(e^x-e^{-x})/(e^x+e^{-x})\), we can obtain \[tanh'(x)=1-tanh^2(x),\] because of \(tanh(x)\in(-1,1)\), we can obtain that \(tanh'(x)\in(0,1]\), and \[|Tanh(x)-Tanh(y)|=|tanh(x)-tanh(y)|\leq 1\cdot |x-y|.\] Therefore, it can be proved that \(Tanh\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=1\).

• Sigmoid(): Given that \(x\), \(y\) and \(Sigmoid(x)=sigmoid(x)=1/(1+e^{-x})\), we can obtain \[sigmoid'(x)=sigmoid(x)\cdot (1-sigmoid(x)),\] because of \(sigmoid(x)\in(0,1)\), we can obtain that \(sigmoid'(x)\in(0,\frac{1}{4}]\), and \[|Sigmoid(x)-Sigmoid(y)|=|sigmoid(x)-sigmoid(y)|\leq \frac{1}{4}\cdot |x-y|.\] Therefore, it can be proved that \(Sigmoid\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=\frac{1}{4}\).

• Softplus(): Given that \(x\), \(y\) and \(Softplus(x)=softplus(x)=\log(1+e^x)\), we can obtain \[softplus'(x)=\frac{e^x}{1+e^x}=\frac{1}{1+e^{-x}}=sigmoid(x),\] because of \(sigmoid(x)\in(0,1)\), we can obtain that \(softplus'(x)\in(0,1)\), and \[|softplus(x)-softplus(y)|=|softplus(x)-softplus(y)|\leq 1\cdot |x-y|.\] Therefore, it can be proved that \(Softplus\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=1\).

• Sum(): Given that \(\mathbf{x}=(x_1, x_2, \dots, x_d)\), \(\mathbf{y}=(y_1, y_2, \dots, y_d)\) and \(Sum(\mathbf{x})=\sum_{i=1}^d x_i\), we then have \[|Sum(\mathbf{x})-Sum(\mathbf{y})|=|\sum_{i=1}^d(x_i-y_i)|\leq \sum_{i=1}^d|x_i-y_i|,\] Based on the Cauchy-Schwarz inequality, for \(\mathbf{a}=(a_1, a_2, \dots, a_d)\) and \(\mathbf{b}=(b_1, b_2, \dots, b_d)\), we have \[\sum_{i=1}^d|a_ib_i|\leq {||a||}_2 {||b||}_2,\] we take \(\mathbf{a}=(1, 1, \dots, 1)\), \(\mathbf{b}=(|x_1-y_1|, |x_2-y_2|, \dots, |x_d-y_d|)\) and we then obtain \[\begin{align} \nonumber \sum_{i=1}^d|x_i-y_i| & =\sum_{i=1}^d|1\cdot (x_i-y_i)|=\sum_{i=1}^d|a_ib_i|\\ & \leq {||a||}_2 {||b||}_2=\sqrt{d}\cdot {||x-y||}_2. \end{align}\] Therefore, it can be proved that \(Sum\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=\sqrt{d}\).

• Mean(): Given that \(\mathbf{x}=(x_1, x_2, \dots, x_d)\), \(\mathbf{y}=(y_1, y_2, \dots, y_d)\) and \(Mean(\mathbf{x})=\sum_{i=1}^d x_i/d\), we then have \[|Mean(\mathbf{x})-Mean(\mathbf{y})|=|\frac{1}{d}\sum_{i=1}^d(x_i-y_i)|\leq \frac{1}{d}\sum_{i=1}^d|x_i-y_i|,\] Based on the Cauchy-Schwarz inequality, for \(\mathbf{a}=(a_1, a_2, \dots, a_d)\) and \(\mathbf{b}=(b_1, b_2, \dots, b_d)\), we have \[\sum_{i=1}^d|a_ib_i|\leq {||a||}_2 {||b||}_2,\] we take \(\mathbf{a}=(1, 1, \dots, 1)\), \(\mathbf{b}=(|x_1-y_1|, |x_2-y_2|, \dots, |x_d-y_d|)\) and we then obtain \[\begin{align} \nonumber \frac{1}{d}\sum_{i=1}^d|x_i-y_i| & =\frac{1}{d}\sum_{i=1}^d|1\cdot (x_i-y_i)|=\frac{1}{d}\sum_{i=1}^d|a_ib_i| \\ & \leq \frac{1}{d}{||a||}_2 {||b||}_2=\frac{1}{\sqrt{d}}\cdot {||x-y||}_2. \end{align}\] Therefore, it can be proved that \(Mean\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=1/\sqrt{d}\).

• FFN(): Given that \(\mathbf{x}\), \(\mathbf{y}\), and \(FFN(\mathbf{x}) = \mathbf{x}\times \mathrm{W}_{ffn}\), we then have \[||FFN(\mathbf{x})-FFN(\mathbf{y})||=||(\mathbf{x}-\mathbf{y})\cdot \mathrm{W}_{ffn}|| \leq ||\mathbf{x}-\mathbf{y}||\cdot ||\mathrm{W}_{ffn}||_2.\] Therefore, it can be proved that \(FFN\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=||\mathrm{W}_{ffn}||_2\).

• FFN_D(): Given that \(\mathbf{x}\), \(\mathbf{y}\), and \(FFN\_D(\mathbf{x})=\mathbf{x}\times \mathrm{W}_{ffnd}\), we then have \[||FFN\_D(\mathbf{x})-FFN\_D(\mathbf{y})||=||(\mathbf{x}-\mathbf{y})\cdot \mathrm{W}_{ffnd}|| \leq ||\mathbf{x}-\mathbf{y}||\cdot ||\mathrm{W}_{ffnd}||_2.\] Therefore, it can be proved that \(FFN\_D\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=||\mathrm{W}_{ffnd}||_2\).

• Add(): Given that \(\mathbf{x}=(x_1,x_2,\dots,x_d)\), \(\mathbf{y}=(y_1,y_2,\dots,y_d)\), \(\mathbf{x'}=(x'_1,x'_2,\dots,x'_d)\), \(\mathbf{y'}=(y'_1,y'_2,\dots,y'_d)\) and \(Add(\mathbf{x},\mathbf{y})=\mathbf{x}+\mathbf{y}\), we then have \[\begin{align} \nonumber & ||Add(\mathbf{x},\mathbf{y})-Add(\mathbf{x}',\mathbf{y}')|| \\ & =||\mathbf{x}-\mathbf{x}'+\mathbf{y}-\mathbf{y}'|| \\ \nonumber & \leq ||\mathbf{x}-\mathbf{x}'|| + ||\mathbf{y}-\mathbf{y}'||. \end{align}\] Therefore, it can be proved that \(Add\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=1\).

• Mul(): Given that \(\mathbf{x}=(x_1,x_2,\dots,x_d)\), \(\mathbf{y}=(y_1,y_2,\dots,y_d)\), \(\mathbf{x'}=(x_1+\delta,x_2+\delta,\dots,x_d+\delta)\), \(\mathbf{y'}=(y_1+\delta,y_2+\delta,\dots,y_d+\delta)\) and \(Mul(\mathbf{x},\mathbf{y})=\mathbf{x}\odot \mathbf{y}\), we then have \[\begin{align} \nonumber & ||Mul(\mathbf{x}, \mathbf{y})-Mul(\mathbf{x}', \mathbf{y}')|| \\ \nonumber & = ||\mathbf{x} \odot \mathbf{y}-\mathbf{x}'\odot \mathbf{y}'|| \\ \nonumber & = ||(x_1y_1,\dots,x_dy_d)-((x_1+\delta)(y_1+\delta),\dots,(x_d+\delta)(y_d+\delta))|| \\ & = ||(-(x_1+y_1)\delta-\delta^2),\dots,(-(x_d+y_d)\delta-\delta^2)|| \\ \nonumber & = ||-\delta((x_1+y_1)+\delta),\dots,-\delta((x_d+y_d)+\delta)||, \end{align}\] because of \(||\mathbf{x}-\mathbf{x}'||+||\mathbf{y}-\mathbf{y}'||=2\sqrt{d}\delta\), we can obtain \[\begin{align} \nonumber & \frac{||Mul(\mathbf{x}, \mathbf{y})-Mul(\mathbf{x}', \mathbf{y}')||}{||\mathbf{x}-\mathbf{x}'||+||\mathbf{y}-\mathbf{y}'||} \\ & = \frac{||-\delta((x_1+y_1)+\delta),\dots,-\delta((x_d+y_d)+\delta)||}{2\sqrt{d}\delta} \\ \nonumber & = \frac{||-((x_1+y_1)+\delta),\dots,-((x_d+y_d)+\delta)||}{2\sqrt{d}}. \end{align}\] Therefore, when \(\mathbf{x},\mathbf{y}\rightarrow \infty\), \(\frac{||Mul(\mathbf{x}, \mathbf{y})-Mul(\mathbf{x}', \mathbf{y}')||}{||\mathbf{x}-\mathbf{x}'||+||\mathbf{y}-\mathbf{y}'||}\rightarrow \infty\). It can be proved that \(Mul\) operator does not satisfy the Lipschitz condition.

• Concat(): Given that \((\mathbf{x}_1,\mathbf{y}_1)\), \((\mathbf{x}_2,\mathbf{y}_2)\), and \(Concat(\mathbf{x},\mathbf{y})=[\mathbf{x} || \mathbf{y}]\times \mathrm{W}_{concat}\), we then have \[\begin{align} \nonumber & ||Concat(\mathbf{x}_1, \mathbf{y}_1)-Concat(\mathbf{x}_2,\mathbf{y}_2)|| \\ \nonumber & =||([\mathbf{x}_1||\mathbf{y}_1]-[\mathbf{x}_2||\mathbf{y}_2]) \cdot \mathrm{W}_{concat}|| \\ & \leq ||([\mathbf{x}_1||\mathbf{y}_1]-[\mathbf{x}_2||\mathbf{y}_2]) \cdot ||\mathrm{W}_{concat}||_2 \\ \nonumber & \leq \sqrt{2} (||\mathbf{x}_1-\mathbf{x}_2||+||\mathbf{y}_1-\mathbf{y}_2||) \cdot ||\mathrm{W}_{concat}||_2. \end{align}\] Therefore, it can be proved that \(Concat\) operator satisfies the Lipschitz condition, and the Lipschitz constant is \(L=\sqrt{2} \cdot ||\mathrm{W}_{concat}||_2\).

7.4 Challenges encountered during supernet training↩︎

We hereby provide a detailed mathematical formalization of the Dimension Incompatibility challenge. Furthermore, regarding our search space design and algorithm, we elucidate the advantages of all-ones padding.

(a) Unary operator mismatch: Operators requiring single input are incorrectly connected to two input nodes. As illustrated in Figure 5 (a), given the triplet \(node_1=\{input_1^L, input_1^R, operator_1\}\), where \(input_1^L=\mathbf{h}_s\in \mathbb{R}^{1\times d}\), \(input_1^R=\mathbf{h}_e\in \mathbb{R}^{1\times d}\), and \(operator_1=Sigmoid(\cdot)\), since unary operators like \(Sigmoid\) only accept a single input, we consistently feed the left input to the operator to ensure proper computation. That is to say, the computation process for the triplet \(node_1\) is \(sigmoid(\mathbf{h_s})\).

(b) Binary operator deficiency: Operators requiring dual inputs are connected to only one input node. As shown in Figure 5 (b), given the triplet \(node_1=\{input_1^L, input_1^R, operator_1\}\), where \(input_1^L=\mathbf{h}_s\in \mathbb{R}^{1\times d}\), \(input_1^R=Zero\), and \(operator_1=Concat(\cdot~,\cdot)\), since \(Zero\) merely serves as a placeholder without actual valid input, while \(Concat\) requires two inputs, we adopt a dimension padding strategy: when the left node provides valid \(d\)-dimensional input, we generate a all-ones padding tensor of identical dimension for the right node. That is to say, the computation process of the triplet \(node_1\) is \(Concat(\mathbf{h}_s,\mathbf{1}^{1\times d})\).

(c) Dimension-sensitive operator left input error: Dimension-sensitive operators receive a scalar as left input. As depicted in Figure 5 (c), the triplet \(node_1=\{input_1^L, input_1^R, operator_1\}\), where \(input_1^L=Sum(\mathbf{h}_s)\in \mathbb{R}^1\), \(input_1^R=\mathbf{h}_e\in \mathbf{R}^{1\times d}\), and \(operator_1=FFN(\cdot)\), has been provided. As previously mentioned, \(FFN\) is a dimension-sensitive operator. When the left input dimension is not \(d\)-dimension, conventional padding methods would conflict with other binary operations. Therefore, we opt to replace the current operator with any non-dimension-sensitive operator to ensure training stability. For example, the computational process of this triplet \(node_1\) then becomes \(Sigmoid(Sum(\mathbf{h}_s))\).

(d) Dimension-sensitive operator right input error: Dimension-sensitive operators receive \(d\)-dimension vectors as left input but the scalar as right input. As can be seen from Figure 5 (d), the left and right input dimensions of the triplet \(node_1=\{input_1^L, input_1^R, operator_1\}\) are inconsistent, where \(input_1^L= Sigmoid(\mathbf{h}_s)\in \mathbb{R}^{1\times d}\), \(input_1^R=Sum(\mathbf{h}_e)\in \mathbb{R}^1\), \(operator_1=Concat(\cdot~,\cdot)\). Given that the \(operator_i\) already meets the scenario (c), we employ an all-ones padding tensor to maintain the original computational path, thereby preventing the omission of critical operations through operator replacement. Other words, \(Concat(Sigmoid(\mathbf{h}_s), [Sum(\mathbf{h}_e)|| \mathbf{1}^{1\times (d-1)}])\) denotes the computation process of the triplet \(node_1\).

It should be emphasized that we consistently employ all-ones padding tensors (rather than zero padding) to address dimension incompatibility issues. All-ones padding tensors minimize interference in computational processes within the defined search space, especially for operations like \(Inv\), \(Mul\), and various activation functions. This leads to better preservation of information flow integrity and gradient stability.

7.5 Implementation details of the sandwich rule↩︎

The Sandwich Rule is an efficient training strategy for weight-sharing NAS frameworks, designed to enhance the quality of supernet optimization. It addresses the performance estimation bias caused by inconsistent optimization difficulty across sub-models of varying sizes. The core methodology involves simultaneously training three types of subnets in each training iteration: the minimal subnet, the maximal subnet and randomly sampled subnets. The minimal subnet, denoted as \(\hat{r}_{ij}=Concat(\mathbf{h}_s,\mathbf{h}_e)\), typically represents the lower bound of the architectural space with the smallest parameter count and simplest computations in the search space. Similarly, the maximal subnet, denoted as \(\hat{r}_{ij}=Concat(Concat(\dots \\(Concat(\mathbf{h}_s,\mathbf{h}_e))),Concat(\dots(Concat(\mathbf{h}_s,\mathbf{h}_c))))\), generally represents the upper bound of the architectural space with the largest parameter count and most complex computations in the search space. The randomly sampled subnets represent stochastic samples from the architectural space.

7.6 Theoretical Analysis of Objective Function↩︎

Theorem A.6. Consider a model with predictive distributions \(\mathcal{P}(y|x;\theta)\) on clean inputs \(x\) and \(\mathcal{P}(y|x+x';\theta)\) on noisy inputs, where the loss function is the cross-entropy loss. Assume that for all \(x\), the two predictive distributions are sufficiently close. Then, there exists a constant \(C>0\) such that the difference between the expected risks under clean and noisy conditions satisfies \[|Risk_{noisy}-Risk_{clean}| \leq C||\mathcal{P}(y|x;\theta)-\mathcal{P}(y|x+x';\theta))||^2.\]

Proof:

For the cross-entropy loss function, we have

\[\ell(Archi_\theta(x),y) = -\log \mathcal{P}(y|x;\theta).\]

Define the KL divergence as

\[D_{KL}(\mathcal{P}(y|x;\theta)||\mathcal{P}(y|x+x');\theta)=\sum_y \mathcal{P}(y|x;\theta)\log \tfrac{\mathcal{P}(y|x;\theta)}{\mathcal{P}(y|x+x';\theta)},\]

and the expected risks under clean and noisy conditions as

\[Risk_{clean}=\mathbb{E}_{x} \mathbb{E}_{y \sim \mathcal{P}(y|x)}[-\log \mathcal{P}(y|x;\theta)],\]

\[Risk_{noisy}=\mathbb{E}_{x} \mathbb{E}_{y \sim \mathcal{P}(y|x+x')}[-\log \mathcal{P}(y|x+x';\theta)].\]

Therefore, we obtain

\[\begin{align} \nonumber & \Delta Risk = Risk_{noisy} - Risk_{clean} \\ \nonumber & =\mathbb{E}_{x} \mathbb{E}_{y \sim \mathcal{P}(y|x+x')}[-\log \mathcal{P}(y|x+x';\theta)] - \mathbb{E}_x \mathbb{E}_{y \sim \mathcal{P}(y|x)}[-\log \mathcal{P}(y|x;\theta)] \\ & =\mathbb{E}_x \sum_y [\mathcal{P}(y|x;\theta)\log \mathcal{P}(y|x;\theta) - \mathcal{P}(y|x+x';\theta)\log \mathcal{P}(y|x+x';\theta)]. \end{align}\]

Then, we introduce the intermediate term \(\sum_y \mathcal{P}(y|x;\theta)\log \mathcal{P}(y|x+x';\theta)\) to decompose the expression

\[\begin{align} \nonumber \Delta Risk \nonumber &= \mathbb{E}_x \sum_y [\mathcal{P}(y|x;\theta)\log \tfrac{\mathcal{P}(y|x;\theta)}{\mathcal{P}(y|x+x';\theta)} \\ \nonumber & \quad + (\mathcal{P}(y|x;\theta)-\mathcal{P}(y|x+x';\theta))\log \mathcal{P}(y|x+x';\theta)] \\ &= \mathbb{E}_x \sum_y [D_{KL}(\mathcal{P}(y|x;\theta)||\mathcal{P}(y|x+x');\theta) \\ \nonumber & \quad + (\mathcal{P}(y|x;\theta)-\mathcal{P}(y|x+x';\theta))\log \mathcal{P}(y|x+x';\theta)]. \end{align}\]

Next, to transform the above equation, we introduce another intermediate term \(\sum_y \mathcal{P}(y|x+x';\theta)\log \mathcal{P}(y|x;\theta)\)

\[\begin{align} \nonumber \Delta Risk & = \mathbb{E}_x \sum_y [D_{KL}(\mathcal{P}(y|x;\theta)||\mathcal{P}(y|x+x');\theta) \\ & \quad - D_{KL}(\mathcal{P}(y|x+x';\theta)||\mathcal{P}(y|x);\theta) \\ \nonumber & \quad + \mathcal{P}(y|x;\theta)\log \mathcal{P}(y|x+x';\theta) - \mathcal{P}(y|x+x';\theta)\log \mathcal{P}(y|x;\theta)]. \end{align}\]

When the two distributions are close, i.e., \(\mathcal{P}(y|x;\theta)\approx \mathcal{P}(y|x+x';\theta)\), we have

\[\mathcal{P}(y|x;\theta)\log \mathcal{P}(y|x+x';\theta) \approx \mathcal{P}(y|x+x';\theta)\log \mathcal{P}(y|x;\theta).\]

Meanwhile, the difference in the bidirectional KL divergence is also a higher-order infinitesimal, \(o(||\mathcal{P}(y|x;\theta)-\mathcal{P}(y|x+x';\theta))||^2)\), when the distributions are close. Therefore, there exists a constant \(C>0\) such that

\[|\Delta Risk| \leq C||\mathcal{P}(y|x;\theta)-\mathcal{P}(y|x+x';\theta))||^2.\]

As the predictive distribution \(\mathcal{P}(y|x;\theta)\) becomes closer to the noisy predictive distribution \(\mathcal{P}(y|x+x';\theta)\), their KL divergence decreases. Since the risk difference \(|\Delta Risk|\) can be upper-bounded by the KL divergence, a smaller KL implies a tighter upper bound on the risk change. This indicates that the model exhibits smaller differences in predictive risk between clean and noisy data, meaning its sensitivity to distributional perturbations is reduced, thereby enhancing robustness.

7.7 Datasets↩︎

We conducted experiments using the following two datasets to evaluate our proposed methods. Table 4 shows the statistics of datasets, and more details of the dataset description are as follows:

• Assistment2009(ASSIST09)¹ is a rich collection of educational data derived from the online learning platform ASSISTments. This dataset spans the academic year 2009-2010 and has been widely used in the CD task.

• SLP-Math² is a publicly available benchmark dataset derived from the Smart Learning Partner (SLP) online platform. SLP collects learning data intentionally from secondary school students across multiple subjects, providing diverse content. SLP-Math specifically focuses on the mathematics subset of these data.

7.8 Baseline Methods↩︎

To validate the effectiveness of our proposed OSCD, we conducted extensive experiments and compared several state-of-the-art baselines. The detailed description of these methods is as follows:

• DINA employs a deterministic inputs, noisy "AND" gate mechanism to establish an association model between learners’ latent abilities and observed outcomes through the mapping relationship between latent attribute vectors and item responses.

• IRT is a widely used cognitive diagnosis method. It employs a linear function to model student ability profiles and exercise characteristics unidimensionally.

• MIRT extends the IRT framework multidimensionally, modeling both learner traits and exercise features across multiple latent dimensions.

• MF is a latent factor factorization model. It decomposes the student-exercise score matrix to obtain latent feature representations, thereby enabling exercise performance prediction.

• NCD innovatively constructs a neural cognitive diagnosis system that maps learners and exercises into a latent space, employing deep neural networks with monotonicity constraints to model their nonlinear interactions, thereby achieving precise assessment of students’ knowledge mastery.

• RCD innovatively employs graph-structured modeling by constructing a triple relational network comprising student-exercise interaction graphs, concept-exercise association graphs, and concept dependency graphs, while designing multi-level attention mechanisms to achieve co-optimization of node features and graph structures, thereby enhancing the predictive performance of cognitive diagnosis.

• KSCD designs a knowledge-sensed cognitive diagnosis approach that maps students, exercises, and interactive/non-interactive concepts into a unified vector space, and generates dynamic mastery vectors through matrix operations, thereby uncovering latent knowledge state association patterns from learning behavior data.

• KaNCD is the extension of NCD. It introduces an explicit knowledge association modeling mechanism that infers mastery levels of uncovered knowledge concepts by mining latent dependency relationships among concepts, significantly enhancing the interpretability of diagnostic results.

• ReliCD incorporates a bayesian uncertainty quantification mechanism that explicitly models knowledge state uncertainty through personalized prior distributions, while optimizing confidence estimation via calibration loss functions, effectively addressing overconfidence issues in diagnosis caused by data noise and sparsity.

• ORCDF explicitly encodes raw response signals as graph-structured features and designs a response-aware graph convolution network to mitigate feature over-smoothing issues in traditional cognitive diagnosis models, significantly improving both prediction accuracy and model interpretability.

• EMO-NAS-CD innovatively reformulates neural architecture search as a multi-objective optimization problem, automatically discovering cognitive diagnosis architectures that achieve both high performance and high interpretability through joint optimization of novel interpretability metrics and model efficiency objectives.

References↩︎

---

1. https://sites.google.com/site/assistmentsdata/home/2009-2010-assistment-data/skill-builder-data-2009-2010↩︎

2. https://aic-fe.bnu.edu.cn/en/data/index.html↩︎