The study employs a secondary analysis of cross-sectional data from the Taiwan Special Needs Education Longitudinal Study (SNELS), initiated by Wang in 2007 and supported by Taiwan's National Science Council [45]. The 2011 SNELS survey was designed in Mandarin Chinese to evaluate the adaptability of students with disabilities in a variety of life and academic contexts. It utilized a nationally representative sample of students with disabilities from across Taiwan, as identified in the Taiwan National Special Education Reporting Network (TNSERN) case list. In SNELS, the sampling was conducted on an individual student basis rather than by schools [45]. This method was chosen to avoid biases in the disability category representation that might result from cluster sampling at the school level. The disability categories for students in SNELS were directly obtained from the TNSERN database, where each student with a disability was diagnosed through school assessments, initial psychiatric evaluation, and a board committee's final diagnosis [46]. The SNELS manual indicates that variables such as the type of disability services, socioeconomic status, school location, and gender were not considered in the sampling framework. This was under the premise that a random sampling process would ensure the sample's distribution of these variables would closely reflect the population's characteristics [45].

Figure 1 illustrates the selection process used for the SNELS 2011 dataset in this analysis, detailing the criteria for inclusion and exclusion. Initially, the study filtered out students with disabilities other than ASD. It then excluded those not in inclusive education settings. Additionally, cases where student data could not be matched with teacher responses were removed. This step is crucial, as SNELS asserts that teachers, who interact directly with the students, offer a more accurate perspective on assessing students' adaptability. The final sample consisted of 508 adolescents (463 boys) with ASD who were enrolled in regular classes within mainstream schools.

Flowchart elucidating the inclusion and exclusion criteria for the participants

In this study, we utilized the SNELS dataset with a focus on the constructs of school bullying victimization (SBV), cumulative risk, and internalizing problems, as outlined in Table 1. We offer a detailed overview of the measures employed for these constructs. To control for potential confounding factors influencing SBV, we included data on the participants' adaptability and social-emotional skills. Recognized as covariates in relation to students' problem behaviors [47, 48], these variables were assessed from the teachers' perspective and served as control variables in our analysis. The reliability of the constructs was determined using ordinal alpha, which is appropriate for the ordinal nature of the item responses, as referenced in [49] and details regarding this are provided later in the methods section. The validity of these constructs is discussed in the results section.

SBV assessment in SNELS includes four items evaluating the extent of bullying endured by participants at school, covering relational, verbal, and physical bullying. Scores were tallied on a four-point Likert scale, with higher scores indicating a greater severity of SBV experienced within the school context. The total of these scores served as the SBV indicator. The ordinal alpha for measuring school bullying victimization was 0.62.

Cumulative risk was assessed via five variables, including the quality of friendship interactions, teacher-student relationships, school connection, experiences of stigma, and the impact of disabilities on learning and daily life. Each variable, except for the quality of friendship interactions (which was assessed with a single item), was measured using multiple items on Likert-type scales, where higher scores indicated greater risk. The scores of all items within each construct were summed to represent that construct. Those scores exceeding the 75th percentile among all participants were identified as risks and coded as ''1''; otherwise, they were coded as ''0''. The five risk indices were combined to form the cumulative risk index (see Table 1). The ordinal alpha of the constructs was 0.66 for teacher-student relationship, 0.73 for school connection, 0.68 for perceived stigma and 0.85 for the impact of disabilities.

Internalizing problems were assessed via seven indicators including mood, sleep disturbances, anxiety, concentration difficulties, reluctance to engage with others, uncontrollable behaviors, and feelings of loneliness and helplessness. All indicators, except mood, were rated on a four-point Likert scale from ''never'' (1) to ''often'' (4), while mood was rated from ''very good'' (1) to ''very bad'' (4). The sum of these item scores served as the overall internalizing problems score, with higher scores signifying greater severity. The ordinal alpha was 0.86 for internalizing problems.

As delineated by Martin et al., adaptability encompasses the cognitive, behavioral, and emotional adjustments individuals make in response to new and uncertain situations [50]. This capacity for adaptation reflects a person's comprehensive ability to navigate the fluctuations of everyday life. A primary aim of the Special Needs Education Longitudinal Study (SNELS) was to assess the adaptability of students with disabilities within their daily routines. Consequently, the survey incorporated items to evaluate the participants' competencies in this area [45]. Four items, detailed in Table 1, were included in the SNELS, prompting teachers to rate the extent of each student's decision-making ability, problem-solving resilience, time management, and task organization on a Likert scale ranging from 1 (strongly agree) to 4 (strongly disagree), with higher scores denoting lower adaptability. The adaptability measure yielded an ordinal alpha of 0.84.

Social-emotional skills, rooted in innate traits and learned from experience, shape how individuals think, feel, and interact, crucial for personal and social development [51]. These skills typically include understanding and managing emotions, solving social problems, and exhibiting positive behaviors [52]. In the SNELS study, teachers assess these skills in students with disabilities using three targeted items (detailed in Table 1), measured on a four-point Likert scale from ''strongly agree'' to ''strongly disagree'', with higher scores indicating lower social-emotional proficiency. The internal consistency for the social-emotional skills, as measured by ordinal alpha, yielded a coefficient of 0.79.

Data analysis was conducted with SPSS 23.0 and the R package within JAMOVI [53]. First, participant demographics were analyzed. Second, confirmatory factor analysis (CFA) was then administered to evaluate the factorial validity and convergent validity of the majority of the study's variables—namely, teacher-student relationships, school connections, perceived stigma, the impact of disabilities, and internalizing problems. The constructs of victimization and the quality of friendship interactions were exempted from this analysis. This is because the validity testing of a single item, such as the quality of friendship interactions, via CFA lacks meaningful interpretation. Moreover, the victimization measure aligns more with causal indicators (i.e., formerly formative measurement), with items typically being independent events, thus rendering them unsuitable for CFA [44, 53].

After assessing the quality of the measurement, we calculated descriptive statistics and bivariate correlations to explore the associations between various risk factors and SBV. It is important to note that the sampling strategy implemented by SNELS was centered on individual students with disabilities, rather than on a cluster sampling approach at the class or school level [45]. Additionally, the intra-class correlation coefficient (ICC) for SBV was found to be negligible, registering close to zero (indeed, it was calculated at -0.14). This low ICC indicates that the variance between two randomly chosen individuals from any class is nearly as significant as the variance between two individuals selected at random from the entire population. These two factors—the predetermined individual-based sampling strategy by SNELS and the minimal ICC—led us to determine that Hierarchical Linear Modeling (HLM) was not the appropriate method for our data analysis. Consequently, we employed hierarchical regression analysis, controlling for variables such as gender, adaptability, and social-emotional skills, to examine the impact of cumulative risk. We also investigated the potential mediating effect of internalizing problems in the relationship between cumulative risk and SBV.

Moreover, in the assessment of Perceived Stigma and the Impact of Disabilities, there were instances of missing data for several response items. To determine the nature of this missing data, Little's Missing Completely at Random (MCAR) test was employed. The results yielded a chi-square value of 27.44 (df = 26) with a p-value of 0.39. This indicates that the data missingness can be considered as completely at random. Given this finding, the Monte Carlo Markov Chain (MCMC) method, available in the PRELIS module of LISREL, was utilized for multiple imputation.