Aflatoxin B₁ (≥ 99%) and Dimethyl sulfoxide (DMSO) were purchased from Shanghai Acmec Biochemical Co., Ltd., and Sangon Biotech Co., Ltd. (Shanghai, China), respectively. Assay kits for the measurements of malondialdehyde (MDA), superoxide dismutase (SOD) Glutathione peroxidase (GSH-Px) and total protein (TP) were purchased from Nanjing Jiancheng Bioengineering Institute (Nanjing, China). Assay kits for the total RNA extraction, PowerUp SYBR Green Master Mix (Thermo Fisher Scientific, America) and HiScript III-RT SuperMix for qPCR (Vazyme, China).

A total of 36 male Wistar rats were purchased from Shanghai SLAC Laboratory Animal Co., Ltd. Body weight (BW): 180-200 g. All rats were housed in an environmentally controlled room at 22-24°C with a relative humidity of 50-55% under a cycle of 12 h each light/dark. Food and water were available ad libitum. All animal procedures were performed in accordance with the Guidelines for Care and Use of Laboratory Animals of Shanghai University of Traditional Chinese Medicine and approved by the Animal Ethics Committee of Animal Center of China and Shanghai University of Traditional Chinese Medicine (Ethics number: PZSHUTCM210115014).

After one week of adaptive feeding, 36 rats were randomly divided into 3 groups (12 in each group): blank group (CK): normal breeding; AFB₁ high-dose (HA) group: AFB₁ was intraperitoneal injection with a dose of 2 mg/kg BW, only once; AFB₁ low-dose (LA) group: AFB₁ was intraperitoneal injection at a dose of 1 mg/kg BW, only once (AFB₁ soluble in 4% DMSO). The experimental doses of AFB₁ were determined according to references and results of pre-experiment (Monmeesil et al., 2019; Deng et al., 2020).

The above dose groups were determined based on references and preliminary experiments. The concentration of DMSO did not affect the normal growth of rats (Chen et al., 2020). Blood samples were collected from ocular venous plexus at the 1st, 4th and 7th day after ASAE and serum were prepared. On the 2nd day after ASAE, 6 rats in each group were randomly selected to weigh their final body weight. The rats were anesthetized with 10% chloral hydrate (350 mg/kg). After anesthesia, the rats were dissected, the blood of abdominal aorta was taken, the liver was separated, washed with normal saline, sucked dry and weighed. Liver tissues of an appropriate size were taken from the same part of the liver of each rat and soaked in 4% paraformaldehyde solution, which was fully fixed for histopathological study, and the rest tissues were frozen in liquid nitrogen and stored in a cryogenic refrigerator at −80°C for use. The remaining 6 rats in each group were treated with the same method on the 7th day after ASAE.

After being fixed in 4% paraformaldehyde for 24 h, the liver was removed and dehydrated in graded alcohol series, transparent in xylene. Tissues were embedded in paraffin, and paraffin blocks were then cut to 4 μm thickness using a microtome (RM2016, Shanghai Leica Instrument Co., Ltd., Germany). These sections were deparaffinized, rehydrated and stained using hematoxylin-eosin (HE) and then analyzed by optical microscopy (Nikon Eclipse E100, Japan) to evaluate histopathological changes. Observe the morphologies of the liver and spleen tissues were observed at 100 ×.

The serum was separated and measured the levels of alanine aminotransferase (ALT), alkaline phosphatase (ALP), aspartate aminotransferase (AST) and total bilirubin (TBIL) by automated chemistry analyzer (ADVIA 2120i, Hitachi, Ltd., Japan). The activities of SOD and GSH-Px, as well as the level of MDA were examined using commercial assay kits.

The apoptosis rate of rat hepatocytes in the model groups were detected by TUNEL method. Paraffin sections were taken in sequence for deparaffinize and rehydrate, antigen retrieval, permeabilization, equilibration at room temperature, TUNEL reaction, BSA blocking, addition of double antibodies, DAPI counterstain in nucleus and then sealed for microscopic examination. Stained nuclei were blue under UV excitation and positive apoptotic cells were green. The number of apoptotic hepatocytes in the total field of view was counted as well as the total number of hepatocytes. Calculation formula: Hepatocyte apoptosis rate = (Number of apoptotic hepatocytes/Total number of hepatocytes) × 100%.

Total RNA in the rats’ liver was extracted using the MagMAX™ mir Vana™ Total RNA Isolation Kit following the specification (Thermo Scientific™ KingFisher™ Flex™, Finland). Total RNA of each sample was quantified and qualified by Agilent 2100/2200 Bioanalyzer (Agilent Technologies, Palo Alto, CA, United States), NanoDrop (Thermo Fisher Scientific Inc.). 1 μg total RNA was used for following library preparation. Illumina RNA sequencing experiments used three parallel samples. The sequencing library construction and Illumina sequencing were conducted at GENEWIZ.

To verify the reliability of the expression profiles observed in the RNA-Seq data, four genes were selected for real-time quantitative PCR (Q-PCR) analysis with the same experimental samples as in the RNA-Seq experiments. Total RNA was extracted using a King Fisher Flex automated nucleic acid extractor and accompanying kit, and reverse transcribed after RNA electrophoresis to determine RNA quality. β-Actin was used as an internal reference and primers designed by Primer 5 are shown in Table 2, and relative expression was calculated according to the 2^–ΔΔct relative quantification formula.

Results were analyzed using GraphPad Prism 8.0 (GraphPad Software, La Jolla, California), SPSS 25 (IBM, NY, United States) and Origin 8.0 software (Origin Inc., Northampton, MA, United States). Analysis of variance (ANOVA) was performed followed by Tukey’s test with a confidence interval of 95% (p ≤ 0.05).

During the experiment, the BW of the rats was weighed before ASAE and on the 1st, 2nd, 4th and 7th after ASAE (Figure 1A). The BW of rats in the CK group showed an upward trend in the whole experimental cycle, while that in the model groups decreased and then increased.

The organ index is one of the main indicators of the biological properties of animals (Zou et al., 2010). Figure 1B shows that there was no significant difference in liver index between rats treated on the 2nd day after ASAE and the CK group. Compared on the 7th day after ASAE with the CK group, the liver coefficient of the model groups was significantly higher (p < 0.05), and the HA group was significantly different from the CK group (p < 0.001).

Representative hepatic histopathological examination results of each group were shown in Figures 1C–G. The rats in the CK group had a normal liver lobular architecture and cell structure. The nucleolus was clear, homogenous cytoplasm (Figure 1C). Livers illustrated different degrees of pathological changes in the model groups. In LA-2d group, the liver exhibited moderate spotty necrosis and degeneration (Figure 1D). In LA-7d group, moderate spotty necrosis, vacuolar degeneration, and mild inflammatory cell infiltration were observed in the portal area (Figure 1E). Hepatocyte lytic necrosis, called hepatocyte bridging necrosis, was widely found in the liver pathological sections of rats in HA-2d group. This phenomenon was common in liver poisoning or severe chronic hepatitis (Figure 1F). In HA-7d group, in addition to severe bridging necrosis of hepatocytes, there were also vacuolar degeneration of hepatocytes and fibrous tissue hyperplasia (Figure 1G).

Transaminase ALT and AST are common indicators to determine whether the liver is damaged. ALT mainly exists in the plasma of hepatocytes, and the elevation of ALT in serum reflects the injury of hepatocytes. AST mainly exists in mitochondria of hepatocytes, and elevated AST indicates more serious liver injury. ALP and TBIL are the detection indexes reflecting hepatic bile metabolism. After ASAE, ALT, AST, ALP and TBIL in the model groups increased first and then decreased (Figure 2). It reached the peak on the 2nd day after ASAE, which was significantly different from that in the CK group (p < 0.001). Except ALP, other indexes returned to normal level on the 7th day after ASAE.

Apoptosis refers to the process by which cells are affected by abnormal factors and controlled by genes to end their life autonomously. As shown in Figures 3A-E, the nuclear membranes of normal hepatocytes were intact and stained blue; the nuclear membranes of apoptotic hepatocytes were clustered or fragmented and stained green. As can be seen from Figure 3F, compared with the CK group, the apoptosis rate in the model groups was significantly increased in a dose-dependent manner. However, the apoptosis rate of hepatocytes on the 7th after ASAE was slightly lower than that on the 2nd after ASAE, and there was still significant difference (p < 0.05).

The content of MDA can reflect the degree of lipid peroxidation and indirectly reflect the degree of cell damage. In this study, the MDA content in model groups were significantly higher than that in CK group (p < 0.05) (Figure 4A), indicating that AFB₁ caused hepatocyte injury. It was also found that MDA content in the liver of the rats after ASAE not only showed a dose-dependent manner, but also first increasing and then decreasing trend with prolonged observation. However, it was still significantly higher than the CK group and did not return to the normal levels on the 7th day (p < 0.05).

Studies have shown that AFB₁ can lead to the imbalance of the body’s oxidative system homeostasis and cause oxidative stress response. The antioxidant enzymes SOD and GSH-Px are the body’s first line of defense against free radicals. In this study, compared with the CK group, the activities of SOD and GSH-Px in the model groups decreased significantly (p < 0.05) (Figures 4B,C). And with the prolonged observation time after ASAE, the activities of SOD and GSH-Px became lower and lower (p < 0.001).

In this study, the DEGs between the model groups and the CK group were obtained by RNA-Seq technology, and a volcano map was drawn (Figures 5A–D). The volcano plot showed that LA-2d, LA-7d, HA-2d and HA-7d groups had 1246, 356, 2354, and 855 DEGs, respectively (FDR ≤ 0.05 and | log₂FC| ≥ 1). From this result, it can be found that the number of DEGs was related to the dose of ASAE and the observation time after ASAE: the higher the dose, the more DEGs; the number of DEGs decreased with the prolonged time of continuous impact. Table 3 shows the classification and number of entries of KEGG pathways enriched in each group of DEGs. Looking at the total number of KEGG pathways in each model group, there were fewer on the 7th day after ASAE than on the 2nd day. The number of entries for each classification in each group shows that the main difference was in the metabolic category. This is related to the fact that the liver is the largest metabolic organ of the body as well as being extremely self-repairing.

In order to explore the genes that were continuously affected and their functions, the screening criteria were strengthened, and the Venn diagram as shown was drawn (FDR ≤ 0.01 and | log₂FC| ≥ 2) (Figure 5E). Compared with CK group under the above filter conditions, the LA-2d, LA-7d, HA-2d and HA-7d groups had 236, 33, 679, and 78 differential genes, respectively. There were 211 common genes in LA-2d and HA-2d, 23 common genes in LA-7d and HA-7d, 20 common genes in LA-2d and LA-7d, and 56 common genes in HA-2d and HA-7d. In order to analyze the accuracy of the results, the first 20 independent DEGs were selected from each group, and all the ones less than 20 were selected for KEGG enrichment analysis. According to KEGG enrichment analysis, these genes are mainly enriched in four types of biological metabolic pathways: Cellular Processes, Environmental Information Processing, Metabolism and Organismal Systems (Tables 4, 5).

By analyzing the DEGs in Tables 4, 5, it was found that DEGs of LA-2d and HA-2d were mainly enriched in pathways related to carbohydrate, lipid and amino acid metabolism, while DEGs of LA-7d and HA-7d were mainly enriched in four pathways, which were related to cell cycle regulation, inflammatory response and oxidation reaction. Lama5, Gtse1 and Fabp4 gene in DEGs were all highly expressed in the model groups. Bcl6 gene was highly expressed on the 7th day after ASAE, and the expression level was higher than that on the 2nd day after ASAE. Lama5 gene enrichment in focal adhesion (ko04510) and metabolic pathways (ko01100). Focal adhesion pathway plays an important role in cell motility, cell proliferation, and cell differentiation (Ye et al., 2020). The dysregulation of focal adhesion is considered to be an essential step in tumor invasion, it can promote tumor invasiveness and metastasis (Shen et al., 2018). Metabolic pathway plays an important role in the metabolism of sugar, fat and protein in the liver. Metabolic dysfunction not only affects the normal growth of the body, but also causes various diseases. Gtse1 gene enrichment in p53 signaling pathway (ko04115). The protein encoded by the p53 gene is a transcription factor that controls the cell cycle. The occurrence of oxidative stress and DNA damage in the body can cause p53 mutation, resulting in cell cycle dysregulation (Lacroix et al., 2020; Bailey et al., 2021). Fabp4 gene enrichment in regulation of lipolysis in adipocytes (ko04923) and PPAR signaling pathways (ko03320). Regulation of lipolysis in adipocytes pathway plays an important role in the dynamic equilibrium of lipid droplet formation and decomposition (Yang and Mottillo, 2020). Dysregulation of this pathway leads to lipid accumulation in the liver causes hepatocyte damage and portal inflammation (An et al., 2021). PPAR signaling pathway plays an essential role in the maintenance of metabolic homeostasis, it can adjust the balance between anabolism and oxidation of adipose tissue to prevent peroxidation (Corrales et al., 2018). Bcl6 gene enrichment in FoxO signaling pathway (ko04068). The pathway can inhibit cell metabolism, growth, differentiation, oxidative stress, and aging, it is suggested to play a pivotal functional role as a tumor suppressor in cancers (Farhan et al., 2017; Lee and Dong, 2017).

Quantitative validation of the predicted differential genes based on transcriptome sequencing data was carried out for Lama5, Ccnd1, Cdkn1a and Gpx2, which were randomly selected. RNA integrity is shown in the plots from electrophoresis experiments using a fully automated nucleic acid electrophoresis analyzer as in Figure 6. The validation results are shown in Table 6. The transcriptome sequenced genes with Q-PCR showed the same trend in different model groups.