LY6D is crucial for lipid accumulation and inflammation in nonalcoholic fatty liver disease

Nonalcoholic fatty liver disease (NAFLD) is a serious metabolic disorder characterized by excess fat accumulation in the liver. Over the past decade, NAFLD prevalence and incidence have risen globally. There are currently no effective licensed drugs for its treatment. Thus, further study is required to identify new targets for NAFLD prevention and treatment. In this study, we fed C57BL6/J mice one of three diets, a standard chow diet, high-sucrose diet, or high-fat diet, and then characterized them. The mice fed a high-sucrose diet had more severely compacted macrovesicular and microvesicular lipid droplets than those in the other groups. Mouse liver transcriptome analysis identified lymphocyte antigen 6 family member D (Ly6d) as a key regulator of hepatic steatosis and the inflammatory response. Data from the Genotype-Tissue Expression project database showed that individuals with high liver Ly6d expression had more severe NAFLD histology than those with low liver Ly6d expression. In AML12 mouse hepatocytes, Ly6d overexpression increased lipid accumulation, while Ly6d knockdown decreased lipid accumulation. Inhibition of Ly6d ameliorated hepatic steatosis in a diet-induced NAFLD mouse model. Western blot analysis showed that Ly6d phosphorylated and activated ATP citrate lyase, which is a key enzyme in de novo lipogenesis. In addition, RNA- and ATAC-sequencing analyses revealed that Ly6d drives NAFLD progression by causing genetic and epigenetic changes. In conclusion, Ly6d is responsible for the regulation of lipid metabolism, and inhibiting Ly6d can prevent diet-induced steatosis in the liver. These findings highlight Ly6d as a novel therapeutic target for NAFLD.


LY6D is crucial for lipid accumulation and inflammation in nonalcoholic fatty
Spearman's rank correlation, since the recorded phenotypes contained non-parametric data 1 .
The R package ggpubr function was applied for the calculation, and the results were visualized using ggplot2. All BXD phenotypes are available on the GeneNetwork website (https://www.genenetwork.org/).

ATAC-seq and data analysis
Single-cell suspensions of AML12 were prepared using trypsin-EDTA, following which the A heatmap of the normalized signal intensity across genomic regions was generated using DeepTools 4 .

Upstream regulator analysis
Upstream regulator analysis was performed using IPA ® , for DEGs from the denoted groups 5 .
The p-value was calculated using Fisher's exact test, with a significance threshold of 0.05 and an activation Z-score threshold of ±2.

Pre-processing for scRNA-seq
The liver tissue was extracted and divided into four pieces, following which the left lobe was stored in fresh DMEM/F12 fluid for an hour. Using pathogen-free scissors, the stored tissues were cut into equal pieces and divided into individual cells. The number of cells were in the range of 0.72×10 5 to 1.88×10 6 . The library used for scRNA-seq was BD WTA (BD Biosciences, San Jose, CA, USA) with 20,000 reads per cell. In the Q30 of FastQC, the quality checking rate was approximately 92% for all samples.

ScRNA-seq
The mitochondrial genome transcript ratio was calculated using the 'PercentageFeatureSet' function. Genes expressed by fewer than three cells were filtered out. Cells with a mitochondrial gene ratio >50% were excluded. Finally, the 'miQC' algorithm 6 was applied to each of the four samples independently, using the linear-type mixture model with a posterior cut-off of 0.75. The normalized data were obtained using the 'NormalizeData' function from the count matrix, while cell cycle scores for the S and G2M phases were estimated using the 'CellCycleScoring' function. We integrated the scRNA-seq data following the procedure suggested by Stuart et al. 7 . Variable features were identified by applying the 'NormalizeData' and 'FindVariableFeatures' functions to the count matrix of each sample independently, from which 2000 integration features were selected using the 'SelectIntegrationFeatures' function.
Anchors for the integration were identified using the 'FindIntegrationAnchors' function with 30 dimensions, to specify the neighbor search space. Finally, integration was performed using the 'IntegrateData' function.
For dimension reduction and analysis following that, the integrated data was scaled using the 'ScaleData' function, after regressing out mitochondrial gene ratio, unique molecular identifier (UMI) counts, S and G2M phase scores. Fifty principal components were calculated using the 'RunPCA' function. The dimension of the principal components was determined as the number of eigenvalues statistically larger than the maximum of those from the permuted matrix obtained by randomly shuffling the scaled expression values for each gene.
For visualization of the integrated data, Uniform Manifold Approximation and Projection (UMAP) was obtained from the principal components space, using the function 'RunUMAP'.
The integrated data were clustered by applying the 'FindNeighbors' and 'FindClusters' functions in series to the principal components space. A resolution value of 0.4 was used to assign the cell type to each cluster using the expression of the marker gene set. In the radar plot, the 'radarchart' function from the 'fmsb' package was used. For the network representation of GO terms, the 'pairwise_termsim' function from the 'enrichplot' package was used to calculate the Jaccard similarity between GO terms and assign it to the associated edge. Edges with a similarity of less than the 90% quantile of all similarity values were eliminated. Finally, the network was imported into Cytoscape (version 3.9.2) 10 , and clustered using the 'GLay' algorithms 11 for community clustering provided by the 'clusterMaker2' plugin 12 .

Generation of loss-of-function Ly6d in mice model, using adeno-associated virus 8 (AAV8)
Seven-week-old C57BL6/J male mice were first habituated to a contingency of 25℃. We then performed Ly6d mRNA knockdown specifically in the mouse liver using an AAV8 virus