The two Chinese rural dogs used in this experiment were purchased from the Experimental Animal Center of Northwest Agriculture and Forest University (Yangling, China). All animal experiments were carried out in accordance with our institute’s laboratory animal requirements and the Guide for the Care and Use of Laboratory Animals (Ministry of Science and Technology of China, 2006).
Cultivation and identification of primary canine BMSCs
The BMSCs were extracted from the long bones of two female 3-month-old Chinese rural dogs, after anesthesia by injecting zoletil (0.1 ml/kg, Virbac Group, France). Followed by the experiment, they were committed to euthanasia with 200 mg/min intravenous infusion of propofol till no heartbeat. The bodies were collected by the Experimental Animal Center of Northwest Agriculture and Forest University (Yangling, China) for harmless disposal. For details about the separation of BMSCs, see our previous paper . The cells were identified by flow cytometry and confirmed to undergo three-lineage differentiation the same as when proliferating in dishes . Specifically, the potential for differentiating into the three lineages was assessed by a 7-day induction procedure. To induce osteogenic differentiation, a total of 2×104 cells were seeded into 12-well plates with α-MEM culture medium consisting of 10% FBS, 100 nmol/L dexamethasone, 30 μg/mL ascorbic acid, and 10 mmol/L β-glycerophosphate (Sigma-Aldrich, St. Louis, MO, USA). The medium for adipogenic differentiation contained 10% FBS, 1 μmol/L dexamethasone, 0.5 mmol/L isobutylmethylxanthine, 1 μg/mL insulin, and 100 μmol/L indomethacin (all from Sigma Aldrich). Chondrogenic differentiation medium was supplied with 10% FBS, 40 ng/mL dexamethasone, 50 μg/mL ascorbic acid, 50 μg/mL L-proline, 1 mmol/L sodium pyruvate (all from Sigma Aldrich), insulin–transferrin–selenium X (Gibco, Carlsbad, CA, USA), and 10 ng/mL transforming growth factor-β3 (PeproTech, Rocky Hill, NJ, USA), which was changed every 2 days.
The induction of BMSCs to transform into IPCs
The protocol used in this study involved three stages, while the BMSCs were placed in a suspended state to form spheroids, which mimicked the characteristics of islets in the pancreas. The dishes used in the protocol were all treated with 2-hydroxyethyl methacrylate (Sigma-Aldrich). In stage 1, passage 4 cells were cultured in Dulbecco’s modified Eagle’s medium (DMEM) high glucose (containing 17.5 mmol/L glucose; Hyclone, Logan, UT, USA), which contained 10 ng/ml basic fibroblast growth factor (bFGF; Invitrogen, Carlsbad, CA, USA), 10 ng/ml epidermal growth factor (EGF; Chemicon, CA, USA), 2% B27 supplement minus insulin (Gibco), 0.5% bovine serum albumin (BSA; Solarbio), and 0.1 mmol/L β-mercaptoethanol. Two days later, the medium was exchanged for stage 2 induction medium and supplied with DMEM high glucose with 10 ng/mL EGF, 20 ng/ml Activin A (Peprotech), 10 mmol/L nicotinamide (Sigma), 2% B27, 0.5% BSA, and 0.1 mmol/L β-mercaptoethanol for 3 days. Stage 3 also lasted 3 days, in which the spheroids were suspended in DMEM low glucose (containing 5.6 mmol/L glucose) with 10 ng/mL EGF, 10 nmol/L exendin-4 (Peprotech), 10 ng/mL betacellulin (Peprotech), 2% B27, 0.5% BSA, and 0.1 mmol/L β-mercaptoethanol. The medium was changed every 2 days.
Collecting those spheroids with a consistent size, we then performed DTZ staining as described previously and GSIS to determine whether a transition occurred after the induction procedure. The selected clusters were precultured in DMEM low glucose for 2 h to remove the insulin present, washed three times with PBS, and incubated in DMEM low glucose for 30 min. The supernatant was then collected and the clusters were transferred into DMEM high glucose stimulated for 30 min. The released insulin was tested by ELISA using Quantikine® ELISA kit, strictly in accordance with the operating instructions.
RNA extraction of cell clusters and quality control
Total RNA was extracted using TRIzol® reagent (Invitrogen, USA) following the manufacturer’s protocol, and RNA degradation and contamination were monitored on 1% agarose gels, as well as the purity being checked using the NanoPhotometer® spectrophotometer (IMPLEN, CA, USA). RNA concentration was measured using Qubit® RNA Assay Kit in a Qubit® 2.0 Fluorometer (Life Technologies, CA, USA). Finally, RNA integrity, assessed by RIN8.0, was determined using the RNA Nano 6000 Assay Kit of the Bioanalyzer 2100 system (Agilent Technologies, CA, USA).
Library construction and sequencing processing
Ribosomal RNA was removed using the Epicentre Ribo-zero™ rRNA Removal Kit (Epicentre, USA), and an rRNA-free residue was obtained by ethanol precipitation. Subsequently, 3 μg of RNA per sample was used for library construction, which was generated using the rRNA-depleted RNA with the NEBNext® Ultra™ Directional RNA Library Prep Kit for Illumina® (NEB, USA), in accordance with the manufacturer’s recommendations. To preferentially select cDNA fragments of 150–200 bp in length, the library fragments were purified with the AMPure XP system (Beckman Coulter, Beverly, USA). PCR was performed with Phusion High-Fidelity DNA polymerase, universal PCR primers, and Index (X) Primer. Finally, products were purified (AMPure XP system) and library quality was assessed on the Agilent Bioanalyzer 2100 system.
Cluster generation was performed on a cBot Cluster Generation System using TruSeq PE Cluster Kit v3-cBot-HS (Illumina), in accordance with the manufacturer’s instructions. Subsequently, the sequencing data were processed on an Illumina Hiseq 4000 platform.
De novo assembly
Clean data were obtained after the processing of raw data through in-house Perl scripts. Q20, Q30, and GC content were calculated, with Q20 and Q30 representing the percentages of bases whose Phred [Phred=−10log10(e)] equaled or exceeded 20 and 30 relative to all bases. For mapping to the reference genome, the genome and gene model annotation files were downloaded directly from NCBI (https://www.ncbi.nlm.nih.gov/genome/?term=dog) and the annotation release ID was 105. The software HISAT2 v2.0.4 was used in this field. Finally, the mapped reads were assembled by StringTie (v1.3.1)  in a reference-based approach.
Identification of DEGs
FPKM was calculated based on the length of the fragments and count of reads mapped to this fragment. This was performed using cuffdiff (v2.1.1) and used to define the expression levels of coding genes. Transcripts with an adjusted P-value <0.05 were considered to be differentially expressed. As thresholds for defining significant differences in gene expression between samples, we used |log2FC| > 1 with a false discovery rate (FDR) < 0.05.
Enrichment analysis of DEGs and PPI analysis
Gene ontology (GO) was applied to annotate the functions of DEGs and lncRNAs, which were divided into three biological modules: molecular function (MF), cellular component (CC), and biological process (BP). This was implemented using the GOseq R package, with correction for bias in gene length . GO terms with a corrected P-value less than 0.05 were considered to be significantly enriched with differentially expressed genes.
Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis was also adopted to define the functions of the DEGs and lncRNAs. We used KOBAS software to test the statistical significance of enrichment of differentially expressed genes and lncRNA target genes in KEGG pathways. The significance formula was the same as that in the GO analysis. The significance, which was calculated using the Q-value derived from multiple hypothesis testing of the P-value, was presented using the magnitude of –log10 (Q-value). Meanwhile, RF (rich factor) was defined as the ratio of DEGs relative to the overall annotated transcripts enriched in the same pathway.
PPI analysis of differentially expressed genes was based on the STRING database, which features known and predicted protein–protein interactions. Canine data were included in the database, so we constructed the networks by extracting the target gene list from the database. We visualized the findings using Cytoscape (v3.6.1).
Lentivirus vectors construction and transfection
The overexpression vector was pCDH-CMV-MCS-CopGFP-T2A-Puro, and complete coding sequences of target genes were synthesized by GeneCreat biological technology Co. (Wuhan, China). The vector used for knockdown was CD513B-U6, and siRNAs were designed by BLOCK-iT™ RNAi Designer (https://rnaidesigner.thermofisher.com/rnaiexpress/). shRNAs derived from siRNAs were synthesized by TsingKe biological technology Co. (Beijing China).
Total RNA was extracted as described above, and cDNA was generated from total RNA. Primers for quantitative real time PCR (qRT-PCR) were designed using Primer Premier 5.0 software (Premier, Canada), and synthesized by TsingKe biological technology Co (Beijing China). qRT-PCR was performed on a Step One Plus™ Real-Time PCR System (Thermo Fisher Scientific, USA). Each 20 μL reaction mixture contained 10 μL of Maxima SYBR Green/ROX qPCR Master Mix (2X) (Thermo Scientific, USA), 0.3 μL of each primer (10 μM), 0.8 μL of cDNA, and 8.6 μL of nuclease-free water. The qRT-PCR run protocol was as follows: 95 °C, 10 min; followed by 40 cycles of 95 °C, 15 s; 60 °C, 30 s; and 72 °C, 15 s in 96-well optical reaction plates. Three biological replicates with three technical replicates for each value determined the Ct values. Expression levels of the tested reference genes were determined by Ct values and calculated by 2-△△Ct.
Bioinformatic analysis of data from the Gene Expression Omnibus (GEO) database using the DAVID database
The gene expression profiles GSE20113 nd GSE52063, chosen from GPL3738, were downloaded from the GEO database; these two series included three normal pancreas samples and four normal bone marrow mesenchymal stem cell samples derived from canine. Because these data were from the same platform, after normalization of the microarray data, we used R (https://www.r-project.org/) to acquire the DEGs, and analyzed these data using the DAVID database (https://david.ncifcrf.gov/).
SPSS 13.0 software was used for analysis of variance and paired T test. Measurement data were presented as mean ± SEM, p < 0.05 was statistically significant difference, expressed with (*); 0.01 < p < 0.05 was expressed with (**); p < 0.01 for the extremely significant difference, expressed with (***); p < 0.001 was expressed with (****).