Proteonomics 2
Alexander Ohlson
10/11/2023
##load package
library(readxl)
library(writexl)
library(limma)## Warning: package 'limma' was built under R version 4.2.2library(ggplot2)
##load data
ms_control_csf <- read_excel("diagnosis_K_csf.xlsx")
head(ms_control_csf)## # A tibble: 6 × 1,011
## ID diagnosis ABHD14B ABL1 ACAN ACP5 ACP6 ACTA2 ACVRL1 ACY1
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 MS_1_K 1 0.520 -0.735 -0.111 -1.58 -1.49 -0.281 -1.56 -1.76
## 2 MS_2_K 1 -1.71 -1.14 0.688 1.20 0.404 -0.858 -0.449 -0.757
## 3 MS_3_K 1 -0.263 -0.312 -0.403 0.0816 -0.989 -0.621 -0.680 -0.329
## 4 MS_4_K 1 -0.468 -0.611 -0.297 -0.293 -1.01 0.162 -0.460 -0.582
## 5 MS_5_K 1 -0.0580 0.251 0.256 1.68 0.473 -1.18 0.381 0.470
## 6 MS_6_K 1 0.481 0.327 -0.191 -0.365 0.369 -0.144 -0.165 -0.174
## # ℹ 1,001 more variables: ADA <dbl>, ADA2 <dbl>, ADAM15 <dbl>, ADAM22 <dbl>,
## # ADAM23 <dbl>, ADAM8 <dbl>, ADAMTS13 <dbl>, ADAMTS15 <dbl>, ADAMTS16 <dbl>,
## # ADAMTS8 <dbl>, ADGRB3 <dbl>, ADGRE2 <dbl>, ADGRE5 <dbl>, ADGRG1 <dbl>,
## # ADGRG2 <dbl>, ADM <dbl>, AGER <dbl>, AGRN <dbl>, AGRP <dbl>, AHCY <dbl>,
## # AKR1B1 <dbl>, AKR1C4 <dbl>, ALCAM <dbl>, ALDH1A1 <dbl>, ALDH3A1 <dbl>,
## # AMBN <dbl>, AMIGO2 <dbl>, AMY2A <dbl>, AMY2B <dbl>, ANG <dbl>,
## # ANGPT1 <dbl>, ANGPT2 <dbl>, ANGPTL1 <dbl>, ANGPTL2 <dbl>, ANGPTL3 <dbl>, …design <- model.matrix(~diagnosis, data = ms_control_csf)
head(design)## (Intercept) diagnosis
## 1 1 1
## 2 1 1
## 3 1 1
## 4 1 1
## 5 1 1
## 6 1 1#transpose the data so that the samples are in different columns
t_ms_control_csf<-as.data.frame(t(ms_control_csf))
#make the name of the samples/ID become your column names
colnames(t_ms_control_csf)<-ms_control_csf$ID
#get rid of the rows that are not protein levels
t_ms_control_csf<-t_ms_control_csf[-c(1,2),]
#make sure all the protein expressions are numeric
t_ms_control_csf[]<-lapply(t_ms_control_csf,as.numeric)
#take a look at your dataframe
head(t_ms_control_csf)## MS_1_K MS_2_K MS_3_K MS_4_K MS_5_K MS_6_K
## ABHD14B 0.5197156 -1.7127358 -0.2628932 -0.4677646 -0.05802174 0.4811136
## ABL1 -0.7347174 -1.1420951 -0.3117541 -0.6113755 0.25094302 0.3267217
## ACAN -0.1106618 0.6878584 -0.4025902 -0.2968317 0.25573450 -0.1907236
## ACP5 -1.5829587 1.1961717 0.0815954 -0.2925337 1.67985571 -0.3650363
## ACP6 -1.4877645 0.4044776 -0.9889641 -1.0088753 0.47314562 0.3686545
## ACTA2 -0.2812268 -0.8580090 -0.6208375 0.1623202 -1.17708912 -0.1439282
## MS_7_K MS_8_K MS_9_K MS_10_K MS_11_K MS_12_K
## ABHD14B -1.2188878 -0.05090515 -0.3717985 0.5788049 -0.27389151 0.6851223
## ABL1 0.8157899 0.83540635 -1.5314686 -0.2854196 -0.37006606 -0.4842715
## ACAN -0.4887703 0.48787874 -0.8572643 1.0374732 0.23825376 -0.5438346
## ACP5 -0.1753680 1.71991236 -0.7495801 0.9716541 1.02492949 -0.3059527
## ACP6 0.5715953 0.31964239 0.1345706 -0.6877438 0.05560663 -1.0902218
## ACTA2 0.5598872 -0.05154042 0.1268194 0.1757935 1.17644754 -0.9623730
## MS_13_K MS_14_K MS_15_K MS_16_K MS_17_K MS_18_K
## ABHD14B 0.05519667 0.26200897 -1.47271274 -0.4968779 0.5581021 -0.02977105
## ABL1 -0.57939798 -1.21948609 -0.49179566 0.3517126 -0.6906475 -1.80609917
## ACAN 1.19270222 -0.85411775 1.67202414 0.1137909 -0.2854692 -1.16142918
## ACP5 0.75414654 -1.84633115 1.88234208 0.7198981 -1.9236405 0.16791748
## ACP6 -0.19711210 -0.73939373 0.19506993 0.1846038 -0.2879887 -0.28569129
## ACTA2 -1.00835306 0.09794819 0.08041161 0.5558239 -1.1454377 -0.09431256
## MS_19_K MS_20_K MS_21_K MS_22_K MS_23_K MS_24_K
## ABHD14B 1.2533709 -0.00302994 0.34180099 0.82184705 -1.03860101 0.4884458
## ABL1 -0.9333006 -1.18992703 0.40814351 -0.99725567 0.78999288 0.8074596
## ACAN 0.2307370 -0.57652358 0.88451676 -0.70255972 0.78312846 2.5712335
## ACP5 -0.4050929 -1.04119249 0.09421324 -1.02116417 0.96103914 0.9736570
## ACP6 -0.7351392 -0.63864660 0.08691993 -0.01450793 -1.07669240 0.6227347
## ACTA2 -1.5703789 -0.25984071 0.20081517 -0.72584311 0.07870073 2.0427971
## MS_25_K MS_26_K MS_27_K MS_28_K MS_29_K MS_30_K
## ABHD14B 0.7004338 -1.4384237 0.5893719 -0.3265111 0.76167956 -0.5952161
## ABL1 0.1824197 -0.6903788 0.6075328 0.8722208 -0.51168085 -1.0010177
## ACAN 0.3476832 -1.0138917 1.7566309 1.0946353 0.53577597 -1.3981184
## ACP5 0.5065964 -0.5753337 0.5532624 1.4747657 -0.32117422 -0.8745568
## ACP6 0.8815799 -1.0438475 1.5051230 0.8395452 0.84218303 -1.2817604
## ACTA2 1.3165263 -0.0429860 0.1625341 -1.1257626 -0.09666502 -0.1464946
## MS_31_K MS_32_K MS_33_K MS_34_K MS_35_K MS_36_K
## ABHD14B 1.1571892 -0.3832281 -0.691829190 0.3851475 0.6859850 -0.3672697
## ABL1 -1.5113147 0.6064580 -0.352061900 -0.7715319 -0.3203531 -2.2271815
## ACAN 0.3639403 0.5817503 -0.009273533 -0.4184977 1.5966821 -1.4933884
## ACP5 0.5081987 0.4102602 -0.391874200 -1.0482024 0.5552653 -0.7325560
## ACP6 -1.6093588 -1.0513354 -0.376397740 0.0675193 0.5923574 -0.1685217
## ACTA2 -0.8605753 -0.5836258 -1.298561985 -0.4741291 -0.5312299 -0.3774641
## MS_37_K MS_38_K MS_39_K MS_40_K MS_41_K MS_42_K
## ABHD14B -0.6312304 -1.4392863 0.60921209 -0.27044104 0.95770915 -0.5921970
## ABL1 -0.5581692 -0.7801309 -0.96500941 -0.43428984 -0.19889220 -1.4105451
## ACAN 0.9558382 -0.6688219 -0.39524830 0.03809928 -0.08216822 0.5775549
## ACP5 -1.2338650 0.6800417 -0.04278050 -0.12569777 -0.27971558 -0.3003448
## ACP6 0.1022362 0.3487433 0.54785506 0.17218059 0.37826973 -0.2989654
## ACTA2 -0.4814004 0.7393163 0.06651067 -0.43349558 -0.11462932 -0.6685285
## MS_43_K MS_44_K MS_45_K MS_46_K MS_47_K MS_48_K
## ABHD14B -0.32909896 -0.2042352 -1.95038661 0.73191929 -0.1089161 0.4880145
## ABL1 -1.05261173 0.4011568 -1.50701518 -0.22952614 -0.4995885 0.2315953
## ACAN 0.45798668 1.4164557 -0.02430697 -0.09755128 0.4385831 -0.5382408
## ACP5 -0.11448190 -1.4111156 -0.32477931 -0.74216960 -0.1587445 -0.7950444
## ACP6 0.87511305 0.4728053 -0.52539117 0.32066348 0.2690986 -0.2494427
## ACTA2 0.09880363 0.4786202 -0.38537693 -0.22583687 -0.7455183 -1.0066422
## MS_49_K MS_50_K MS_51_K CONTROL_1_K CONTROL_2_K CONTROL_3_K
## ABHD14B -2.0972471 -0.7466053 -0.3032205 0.1600046 -0.1095631 0.41124162
## ABL1 -0.6164812 -0.3582424 -0.6992466 0.7056152 -0.2085661 0.71260187
## ACAN -1.4364012 0.5207425 -0.4018910 0.9012983 -1.5960004 0.17235137
## ACP5 -0.9190197 -1.0864565 -0.6065778 0.3315489 -0.9444557 -0.12309408
## ACP6 -0.4685507 -0.5466638 -0.1031722 0.3610815 -1.1293634 0.04352378
## ACTA2 -0.7463737 0.6687423 -0.9587374 -0.5406398 -0.8823891 -0.07698984
## CONTROL_4_K CONTROL_5_K CONTROL_6_K CONTROL_7_K CONTROL_8_K CONTROL_9_K
## ABHD14B 1.9089595 -1.3211079 -0.1964717 0.1183833 -1.3114035 0.43453226
## ABL1 0.1450678 1.7237906 0.1477550 0.5535204 -0.7718006 -1.37829888
## ACAN 0.8247326 1.0484861 0.4745934 -0.5987241 -0.5901586 0.39488120
## ACP5 0.7473369 0.7299123 1.2139969 -0.6854894 -0.5623153 -0.04778758
## ACP6 0.3434677 0.1429095 1.1639952 0.8550317 0.2208523 0.64987855
## ACTA2 0.2591991 0.2224151 0.1139877 0.4576618 0.5889723 0.50535276
## CONTROL_10_K CONTROL_11_K CONTROL_12_K CONTROL_13_K CONTROL_14_K
## ABHD14B 1.15417005 -0.8634899 0.75607255 -1.56716924 1.075456
## ABL1 -0.99886798 1.0552183 -0.01293881 -1.14908177 2.294818
## ACAN -0.04703193 0.9942958 -0.90236459 -1.64372282 2.328426
## ACP5 -0.24166177 0.2744682 -0.28672550 0.06156707 2.609971
## ACP6 -1.57540770 1.1014537 -0.05126701 -0.29394507 1.675900
## ACTA2 -0.32656524 0.9873947 -1.43094174 -2.22308166 2.420475
## CONTROL_15_K CONTROL_16_K CONTROL_17_K CONTROL_18_K CONTROL_19_K
## ABHD14B -1.234846250 -0.4722933 -0.0403381 -0.43994518 1.7825861
## ABL1 0.015007940 -1.4059769 -0.2464554 0.34768178 -0.5310286
## ACAN 0.468300312 -0.7357731 -0.7079788 0.03757485 0.3670868
## ACP5 1.186357790 -0.9060013 -0.5679232 -0.94405512 0.1845410
## ACP6 -0.265609940 -1.5916600 -0.9183390 -0.27258735 1.2992040
## ACTA2 0.008340566 -1.3770489 0.8774703 0.09901749 0.5620258
## CONTROL_20_K
## ABHD14B -0.2447782
## ABL1 -0.9980618
## ACAN -1.1902724
## ACP5 -2.0301912
## ACP6 -1.7591181
## ACTA2 -0.2446566fit <- lmFit(t_ms_control_csf, design)
ebayes <- eBayes(fit)
results <- topTable(ebayes, coef = 2, number = Inf, adjust = "BH")
results$proteins<-rownames(results)
head(results)## logFC AveExpr t P.Value adj.P.Val B
## MZB1 3.140727 2.256015 10.225869 1.015850e-16 1.024993e-13 27.554535
## TNFRSF13B 2.564393 1.842028 8.454328 4.819481e-13 2.431428e-10 19.336957
## IL12B 4.902686 3.521648 8.225389 1.433047e-12 4.819815e-10 18.278657
## CD79B 2.694839 1.935729 8.101793 2.577291e-12 6.501215e-10 17.708697
## CD27 2.027617 1.456457 8.050987 3.279493e-12 6.618018e-10 17.474730
## IL12A_IL12B 2.435523 1.749460 6.228306 1.495415e-08 2.514789e-06 9.308632
## proteins
## MZB1 MZB1
## TNFRSF13B TNFRSF13B
## IL12B IL12B
## CD79B CD79B
## CD27 CD27
## IL12A_IL12B IL12A_IL12Bnrow(results)## [1] 1009significant_results <- results[results$adj.P.Val < 0.05, ]
significant_results$proteins<-rownames(significant_results)
head(significant_results)## logFC AveExpr t P.Value adj.P.Val B
## MZB1 3.140727 2.256015 10.225869 1.015850e-16 1.024993e-13 27.554535
## TNFRSF13B 2.564393 1.842028 8.454328 4.819481e-13 2.431428e-10 19.336957
## IL12B 4.902686 3.521648 8.225389 1.433047e-12 4.819815e-10 18.278657
## CD79B 2.694839 1.935729 8.101793 2.577291e-12 6.501215e-10 17.708697
## CD27 2.027617 1.456457 8.050987 3.279493e-12 6.618018e-10 17.474730
## IL12A_IL12B 2.435523 1.749460 6.228306 1.495415e-08 2.514789e-06 9.308632
## proteins
## MZB1 MZB1
## TNFRSF13B TNFRSF13B
## IL12B IL12B
## CD79B CD79B
## CD27 CD27
## IL12A_IL12B IL12A_IL12Bggplot(data = results, aes(x = logFC, y = -log10(P.Value), col = ifelse(logFC > 0, "Upregulated", "Downregulated"))) +
labs(color = "Regulation")+geom_point(position = position_nudge(x = 0.0),size= 2) +
theme(axis.title.x=element_text(size=20),axis.text.x=element_text(size=20),axis.text.y=element_text(size=20),axis.title.y=element_text(size=20),axis.line = element_line(color = "black", size =0.5),legend.title=element_text(size=20),legend.text = element_text(size = 20),panel.grid.major = element_blank(),panel.grid.minor = element_blank(), panel.background = element_blank()) + geom_hline(yintercept = -log10(0.05), linetype = "solid", color = "black") +
geom_vline(xintercept = 0, linetype = "solid", color = "black")+
geom_text(data = subset(results, adj.P.Val < 0.000001), aes(label = proteins), nudge_x = 0.5, nudge_y = 0.5,size=6) +
scale_color_manual(values = c("Upregulated" = "red", "Downregulated" = "blue")) +
geom_point(data = subset(results, P.Value > 0.05), aes(col = "gray50"), size = 2)+
coord_cartesian(xlim = c(-1, 7)) +
scale_x_continuous(breaks = seq(-1, 7, by = 1), labels = seq(-1, 7, by = 1))+guides(color = guide_legend(override.aes = list(size = 5)))## Warning: The `size` argument of `element_line()` is deprecated as of ggplot2 3.4.0.
## ℹ Please use the `linewidth` argument instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
write_xlsx(significant_results,"/Volumes/KODEN/2023/significant_results_DEA_MS_control_L.xlsx")
knitr::opts_chunk$set(echo = TRUE)R Markdown
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:
summary(cars)## speed dist
## Min. : 4.0 Min. : 2.00
## 1st Qu.:12.0 1st Qu.: 26.00
## Median :15.0 Median : 36.00
## Mean :15.4 Mean : 42.98
## 3rd Qu.:19.0 3rd Qu.: 56.00
## Max. :25.0 Max. :120.00Including Plots
You can also embed plots, for example:
Note that the echo = FALSE parameter was added to the
code chunk to prevent printing of the R code that generated the
plot.