setwd("/home/mchopra/Documents/PhD-Year1/PRS/")
library(RegParallel)
## Loading required package: doParallel
## Loading required package: foreach
## Loading required package: iterators
## Loading required package: parallel
## Loading required package: data.table
## Loading required package: stringr
## Loading required package: survival
## Loading required package: arm
## Loading required package: MASS
## Loading required package: Matrix
## Loading required package: lme4
##
## arm (Version 1.13-1, built: 2022-8-25)
## Working directory is /home/mchopra/Documents/PhD-Year1/PRS
library(ggplot2)
GWAS of cellular proportion in GTEx samples and calculating PGS on UKBB data. Further, need to find if cellular proportion is associated with the distensibility. Can perform correlations as well.
aao_dist <- read.table("/home/mchopra/Documents/PhD-Year1/GWAS/norm_pheno_aao.txt", header = T, sep = ",")
dao_dist <- read.table("/home/mchopra/Documents/PhD-Year1/GWAS/norm_pheno_dao.txt", header = T, sep = ",")
rownames(aao_dist) <- aao_dist[,1]
rownames(dao_dist) <- dao_dist[,1]
intersect <- intersect(rownames(aao_dist), rownames(dao_dist))
aao_dist <- aao_dist[intersect, ]
dao_dist <- dao_dist[intersect,]
distensibility <- merge(aao_dist, dao_dist, by = "IID", all = TRUE)
#for VSMCI
##Case1: where PGS was calculated (without additional optimisation)
case1 <- read.table("cellp_pgs/VSMCI_nopheno_nocovar.all_score", header = T)
##distensibility_results:
#distensibility
#i want both of them in one table
#making rows
rownames(case1) <- case1[,1]
rownames(distensibility) <- distensibility[,1]
intersect <- intersect(rownames(case1), rownames(distensibility))
case1 <- case1[intersect, ]
distensibility <- distensibility[intersect,]
case1_data <- merge(case1, distensibility, by = "IID", all = TRUE)
colnames(case1_data)
## [1] "IID" "FID" "Pt_0.001"
## [4] "Pt_0.05" "Pt_0.1" "Pt_0.2"
## [7] "Pt_0.3" "Pt_0.4" "Pt_0.5"
## [10] "Pt_1" "FID.x" "AAo_distensibility"
## [13] "FID.y" "DAo_distensibility"
# Define the response variables
distens <- c("AAo_distensibility", "DAo_distensibility")
# Initialize a list to store the model results
macro_res_list <- list()
# Use a for loop to fit linear regression models for each response variable
for (ctypes in distens) {
formula <- paste(ctypes, '~ as.numeric([*])')
# Note: Replace 'as.numeric([*])' with the actual predictor variable
# you want to use. I've used a placeholder here.
macro_res <- RegParallel(
data = case1_data,
formula = formula,
FUN = function(formula, data)
lm(formula = formula, data = data),
FUNtype = 'lm',
variables = colnames(case1_data)[3:10],
blocksize = 8,
p.adjust = "BH",
cores = 2,
nestedParallel = FALSE,
conflevel = 95
)
# Store the results in the list
macro_res_list[[ctypes]] <- macro_res
}
##
## ##############################
## #RegParallel
## ##############################
## System is:
## -- Linux
## Blocksize:
## -- 8
## Cores / Threads:
## -- 2
## Terms included in model:
## -- AAo_distensibility
## First 5 formulae:
## -- AAo_distensibility ~ as.numeric(Pt_0.001)
## -- AAo_distensibility ~ as.numeric(Pt_0.05)
## -- AAo_distensibility ~ as.numeric(Pt_0.1)
## -- AAo_distensibility ~ as.numeric(Pt_0.2)
## -- AAo_distensibility ~ as.numeric(Pt_0.3)
## Processing 8 formulae, batch 1 of 1
## -- index1: 1; index2: 8
## Done!
##
## ##############################
## #RegParallel
## ##############################
## System is:
## -- Linux
## Blocksize:
## -- 8
## Cores / Threads:
## -- 2
## Terms included in model:
## -- DAo_distensibility
## First 5 formulae:
## -- DAo_distensibility ~ as.numeric(Pt_0.001)
## -- DAo_distensibility ~ as.numeric(Pt_0.05)
## -- DAo_distensibility ~ as.numeric(Pt_0.1)
## -- DAo_distensibility ~ as.numeric(Pt_0.2)
## -- DAo_distensibility ~ as.numeric(Pt_0.3)
## Processing 8 formulae, batch 1 of 1
## -- index1: 1; index2: 8
## Done!
PGS_results <- dplyr::bind_rows(macro_res_list, .id ="distens")
PGS_results
## distens Variable Term Beta StandardError
## <char> <char> <char> <num> <num>
## 1: AAo_distensibility Pt_0.001 as.numeric(Pt_0.001) 39.57533 25.36870
## 2: AAo_distensibility Pt_0.05 as.numeric(Pt_0.05) 455.94933 221.70091
## 3: AAo_distensibility Pt_0.1 as.numeric(Pt_0.1) 414.04947 323.68739
## 4: AAo_distensibility Pt_0.2 as.numeric(Pt_0.2) 297.46584 475.53976
## 5: AAo_distensibility Pt_0.3 as.numeric(Pt_0.3) 173.79724 615.44473
## 6: AAo_distensibility Pt_0.4 as.numeric(Pt_0.4) 136.69057 728.77421
## 7: AAo_distensibility Pt_0.5 as.numeric(Pt_0.5) 157.12139 843.76856
## 8: AAo_distensibility Pt_1 as.numeric(Pt_1) 230.04442 1315.66968
## 9: DAo_distensibility Pt_0.001 as.numeric(Pt_0.001) 16.07625 25.67877
## 10: DAo_distensibility Pt_0.05 as.numeric(Pt_0.05) 244.08982 224.41310
## 11: DAo_distensibility Pt_0.1 as.numeric(Pt_0.1) 39.11156 327.64227
## 12: DAo_distensibility Pt_0.2 as.numeric(Pt_0.2) -397.13080 481.34023
## 13: DAo_distensibility Pt_0.3 as.numeric(Pt_0.3) -504.98462 622.94981
## 14: DAo_distensibility Pt_0.4 as.numeric(Pt_0.4) -747.25118 737.65805
## 15: DAo_distensibility Pt_0.5 as.numeric(Pt_0.5) -605.49052 854.05894
## 16: DAo_distensibility Pt_1 as.numeric(Pt_1) -1007.42858 1331.71419
## t P OR ORlower ORupper P.adjust
## <num> <num> <num> <num> <num> <num>
## 1: 1.5600062 0.11876514 1.539379e+17 3.921672e-05 6.042544e+38 0.4750605
## 2: 2.0565966 0.03973055 1.038191e+198 2.014286e+09 NA 0.3178444
## 3: 1.2791647 0.20084539 6.597791e+179 1.977857e-96 NA 0.5355877
## 4: 0.6255331 0.53162429 1.540896e+129 0.000000e+00 NA 0.8611985
## 5: 0.2823929 0.77764348 3.014262e+75 0.000000e+00 NA 0.8611985
## 6: 0.1875623 0.85122059 2.311858e+59 0.000000e+00 NA 0.8611985
## 7: 0.1862138 0.85227787 1.725644e+68 0.000000e+00 NA 0.8611985
## 8: 0.1748497 0.86119850 8.072764e+99 0.000000e+00 NA 0.8611985
## 9: 0.6260522 0.53128375 9.590191e+06 1.330506e-15 6.912542e+28 0.6071814
## 10: 1.0876808 0.27674162 1.015921e+106 9.684669e-86 1.065701e+297 0.6071814
## 11: 0.1193727 0.90498057 9.681282e+16 1.248319e-262 7.508276e+295 0.9049806
## 12: -0.8250522 0.40934616 3.375093e-173 0.000000e+00 NA 0.6071814
## 13: -0.8106345 0.41757968 4.874880e-220 0.000000e+00 NA 0.6071814
## 14: -1.0130049 0.31106305 0.000000e+00 0.000000e+00 NA 0.6071814
## 15: -0.7089564 0.47835507 1.093476e-263 0.000000e+00 NA 0.6071814
## 16: -0.7564901 0.44935921 0.000000e+00 0.000000e+00 NA 0.6071814
# List of variables to fit the linear model
variables <- c("Pt_0.001","Pt_0.05","Pt_0.1","Pt_0.2","Pt_0.3","Pt_0.4","Pt_0.5","Pt_1")
# Empty list to store the model summaries
model_summaries <- list()
# Loop through each variable and fit a linear model
for (var in variables) {
# Construct the formula dynamically
formula <- as.formula(paste("AAo_distensibility ~", var))
# Fit the linear model
fit <- lm(formula, data = case1_data)
# Store the summary of the model
model_summary <- summary(fit)
model_summaries[[var]] <- model_summary
}
# Access the summaries of each model
for (var in variables) {
print(var)
print(model_summaries[[var]])
}
## [1] "Pt_0.001"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.9992 -0.6616 0.0042 0.6695 4.3157
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.039440 0.006904 -5.712 1.12e-08 ***
## Pt_0.001 39.575326 25.368697 1.560 0.119
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 5.174e-05, Adjusted R-squared: 3.048e-05
## F-statistic: 2.434 on 1 and 47037 DF, p-value: 0.1188
##
## [1] "Pt_0.05"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0013 -0.6617 0.0041 0.6717 4.3344
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.023e-03 1.966e-02 0.408 0.6832
## Pt_0.05 4.559e+02 2.217e+02 2.057 0.0397 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 8.991e-05, Adjusted R-squared: 6.865e-05
## F-statistic: 4.23 on 1 and 47037 DF, p-value: 0.03973
##
## [1] "Pt_0.1"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0065 -0.6620 0.0037 0.6710 4.3239
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.009703 0.017500 -0.554 0.579
## Pt_0.1 414.049467 323.687387 1.279 0.201
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 3.479e-05, Adjusted R-squared: 1.353e-05
## F-statistic: 1.636 on 1 and 47037 DF, p-value: 0.2008
##
## [1] "Pt_0.2"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0080 -0.6621 0.0042 0.6706 4.3240
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.01408 0.02794 -0.504 0.614
## Pt_0.2 297.46584 475.53976 0.626 0.532
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 8.319e-06, Adjusted R-squared: -1.294e-05
## F-statistic: 0.3913 on 1 and 47037 DF, p-value: 0.5316
##
## [1] "Pt_0.3"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0079 -0.6622 0.0045 0.6705 4.3233
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.02370 0.02737 -0.866 0.387
## Pt_0.3 173.79724 615.44473 0.282 0.778
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 1.695e-06, Adjusted R-squared: -1.956e-05
## F-statistic: 0.07975 on 1 and 47037 DF, p-value: 0.7776
##
## [1] "Pt_0.4"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0087 -0.6620 0.0043 0.6703 4.3228
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.02612 0.02811 -0.929 0.353
## Pt_0.4 136.69057 728.77421 0.188 0.851
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 7.479e-07, Adjusted R-squared: -2.051e-05
## F-statistic: 0.03518 on 1 and 47037 DF, p-value: 0.8512
##
## [1] "Pt_0.5"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0089 -0.6620 0.0043 0.6703 4.3227
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.02608 0.02853 -0.914 0.361
## Pt_0.5 157.12139 843.76856 0.186 0.852
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 7.372e-07, Adjusted R-squared: -2.052e-05
## F-statistic: 0.03468 on 1 and 47037 DF, p-value: 0.8523
##
## [1] "Pt_1"
##
## Call:
## lm(formula = formula, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0089 -0.6620 0.0043 0.6703 4.3227
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.02650 0.02794 -0.948 0.343
## Pt_1 230.04442 1315.66968 0.175 0.861
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 6.5e-07, Adjusted R-squared: -2.061e-05
## F-statistic: 0.03057 on 1 and 47037 DF, p-value: 0.8612
the below shows that AAo_distensibility has relationship with VSMCI
fit <- lm(AAo_distensibility ~ Pt_0.05, data = case1_data)
summary(fit)
##
## Call:
## lm(formula = AAo_distensibility ~ Pt_0.05, data = case1_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0013 -0.6617 0.0041 0.6717 4.3344
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.023e-03 1.966e-02 0.408 0.6832
## Pt_0.05 4.559e+02 2.217e+02 2.057 0.0397 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9841 on 47037 degrees of freedom
## Multiple R-squared: 8.991e-05, Adjusted R-squared: 6.865e-05
## F-statistic: 4.23 on 1 and 47037 DF, p-value: 0.03973
ggplot(case1_data, aes(x = AAo_distensibility, y = Pt_0.05)) +
geom_point(color = "lightblue") +
geom_smooth(method = "lm", color = "blue")
## `geom_smooth()` using formula = 'y ~ x'
OBSERVATION:: PRS of cellular proportion (VSMCI) shows some association
WITH AAo_distensibility
CORRELATION:
library(corrplot)
## corrplot 0.92 loaded
##
## Attaching package: 'corrplot'
## The following object is masked from 'package:arm':
##
## corrplot
x <- case1_data %>% dplyr::select(AAo_distensibility, DAo_distensibility, Pt_0.001, Pt_0.05, Pt_0.1, Pt_0.2, Pt_0.3, Pt_0.4, Pt_0.5, Pt_1)
# Calculate correlation matrix
x <- cor(x)
x
## AAo_distensibility DAo_distensibility Pt_0.001
## AAo_distensibility 1.0000000000 0.7081264374 0.007192755
## DAo_distensibility 0.7081264374 1.0000000000 0.002886615
## Pt_0.001 0.0071927554 0.0028866154 1.000000000
## Pt_0.05 0.0094822145 0.0050150606 0.294050659
## Pt_0.1 0.0058979229 0.0005504087 0.257896076
## Pt_0.2 0.0028842218 -0.0038041572 0.234180869
## Pt_0.3 0.0013020679 -0.0037376812 0.223085204
## Pt_0.4 0.0008648197 -0.0046707539 0.215927973
## Pt_0.5 0.0008586022 -0.0032688679 0.213198984
## Pt_1 0.0008062039 -0.0034880346 0.210561907
## Pt_0.05 Pt_0.1 Pt_0.2 Pt_0.3
## AAo_distensibility 0.009482214 0.0058979229 0.002884222 0.001302068
## DAo_distensibility 0.005015061 0.0005504087 -0.003804157 -0.003737681
## Pt_0.001 0.294050659 0.2578960755 0.234180869 0.223085204
## Pt_0.05 1.000000000 0.8812401402 0.801823167 0.764512573
## Pt_0.1 0.881240140 1.0000000000 0.911640737 0.870918970
## Pt_0.2 0.801823167 0.9116407375 1.000000000 0.957840848
## Pt_0.3 0.764512573 0.8709189704 0.957840848 1.000000000
## Pt_0.4 0.745580956 0.8510620822 0.938394868 0.980802448
## Pt_0.5 0.735226314 0.8396758219 0.926360802 0.969175951
## Pt_1 0.724614189 0.8287193460 0.914734182 0.957287731
## Pt_0.4 Pt_0.5 Pt_1
## AAo_distensibility 0.0008648197 0.0008586022 0.0008062039
## DAo_distensibility -0.0046707539 -0.0032688679 -0.0034880346
## Pt_0.001 0.2159279729 0.2131989836 0.2105619068
## Pt_0.05 0.7455809563 0.7352263138 0.7246141889
## Pt_0.1 0.8510620822 0.8396758219 0.8287193460
## Pt_0.2 0.9383948679 0.9263608017 0.9147341820
## Pt_0.3 0.9808024476 0.9691759512 0.9572877312
## Pt_0.4 1.0000000000 0.9887484918 0.9772486191
## Pt_0.5 0.9887484918 1.0000000000 0.9885711126
## Pt_1 0.9772486191 0.9885711126 1.0000000000
# Create correlation plot
corrplot(x, method = 'number', addCoef.col = "black", tl.col = "black")
corrplot(x, method="color",
type="upper", order="hclust",
addCoef.col = "black", # Add coefficient of correlation
tl.col="black", tl.srt=45, #Text label color and rotation
# hide correlation coefficient on the principal diagonal
diag=FALSE
)
…. NO CORRELATION FROM THE ABOVE ANALYSIS. ….
Now, to the other way around where we are going to do the same but the pgs is calculated for distensibility in GTEx data. Therefore, the base data is UKB and target data is GTEx.
#PGS of distensibility
pgs_aao <- read.table("dis_pgs/aao.all_score", header = T)
#cell proportion file
cell_prop <- read.csv("/home/mchopra/Documents/PhD-Year1/deconvolution/Deconvolution_results/results/CIBERSORTx_Results.csv")
#Round1
rownames(pgs_aao) <- pgs_aao[,1]
rownames(cell_prop) <- cell_prop[,1]
intersect <- intersect(rownames(cell_prop), rownames(pgs_aao))
cell_prop <- cell_prop[intersect, ]
pgs_aao <- pgs_aao[intersect,]
#changing the column1 name
names(pgs_aao)[1] <- "Mixture"
case2_aao <- merge(cell_prop, pgs_aao, by = "Mixture", all = TRUE)
# List of variables to fit the linear model
variables <- c("Pt_0.001","Pt_0.05","Pt_0.1","Pt_0.2","Pt_0.3","Pt_0.4","Pt_0.5","Pt_1")
# Empty list to store the model summaries
model_summaries <- list()
# Loop through each variable and fit a linear model
for (var in variables) {
# Construct the formula dynamically
formula <- as.formula(paste("VSMC_I ~", var))
# Fit the linear model
fit <- lm(formula, data = case2_aao)
# Store the summary of the model
model_summary <- summary(fit)
model_summaries[[var]] <- model_summary
}
# Access the summaries of each model
for (var in variables) {
print(var)
print(model_summaries[[var]])
}
## [1] "Pt_0.001"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.27423 -0.06804 0.00040 0.07529 0.35457
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.25973 0.00767 33.862 <2e-16 ***
## Pt_0.001 -23.00487 29.02425 -0.793 0.428
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.113 on 385 degrees of freedom
## Multiple R-squared: 0.001629, Adjusted R-squared: -0.0009641
## F-statistic: 0.6282 on 1 and 385 DF, p-value: 0.4285
##
## [1] "Pt_0.05"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.27594 -0.06560 0.00131 0.07568 0.35608
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.31097 0.03044 10.216 <2e-16 ***
## Pt_0.05 -329.17730 208.45690 -1.579 0.115
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1127 on 385 degrees of freedom
## Multiple R-squared: 0.006435, Adjusted R-squared: 0.003855
## F-statistic: 2.494 on 1 and 385 DF, p-value: 0.1151
##
## [1] "Pt_0.1"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.26901 -0.06960 -0.00138 0.07530 0.35532
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.2889 0.0236 12.244 <2e-16 ***
## Pt_0.1 -335.6677 305.6644 -1.098 0.273
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1129 on 385 degrees of freedom
## Multiple R-squared: 0.003123, Adjusted R-squared: 0.0005333
## F-statistic: 1.206 on 1 and 385 DF, p-value: 0.2728
##
## [1] "Pt_0.2"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.27049 -0.06828 -0.00275 0.07488 0.35697
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.28046 0.01471 19.066 <2e-16 ***
## Pt_0.2 -494.00773 400.82453 -1.232 0.219
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1129 on 385 degrees of freedom
## Multiple R-squared: 0.00393, Adjusted R-squared: 0.001343
## F-statistic: 1.519 on 1 and 385 DF, p-value: 0.2185
##
## [1] "Pt_0.3"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.27475 -0.06914 -0.00143 0.07415 0.35713
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.27838 0.01217 22.884 <2e-16 ***
## Pt_0.3 -649.61454 476.74683 -1.363 0.174
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1128 on 385 degrees of freedom
## Multiple R-squared: 0.004799, Adjusted R-squared: 0.002214
## F-statistic: 1.857 on 1 and 385 DF, p-value: 0.1738
##
## [1] "Pt_0.4"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.27501 -0.06763 -0.00154 0.07356 0.35804
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.27493 0.01011 27.18 <2e-16 ***
## Pt_0.4 -737.08929 550.20681 -1.34 0.181
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1128 on 385 degrees of freedom
## Multiple R-squared: 0.00464, Adjusted R-squared: 0.002055
## F-statistic: 1.795 on 1 and 385 DF, p-value: 0.1811
##
## [1] "Pt_0.5"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.27476 -0.06786 -0.00198 0.07303 0.35830
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.724e-01 8.792e-03 30.981 <2e-16 ***
## Pt_0.5 -7.912e+02 6.117e+02 -1.293 0.197
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1128 on 385 degrees of freedom
## Multiple R-squared: 0.004326, Adjusted R-squared: 0.00174
## F-statistic: 1.673 on 1 and 385 DF, p-value: 0.1967
##
## [1] "Pt_1"
##
## Call:
## lm(formula = formula, data = case2_aao)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.27633 -0.06791 -0.00101 0.07379 0.35902
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.734e-01 8.566e-03 31.916 <2e-16 ***
## Pt_1 -1.228e+03 8.113e+02 -1.513 0.131
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1127 on 385 degrees of freedom
## Multiple R-squared: 0.005914, Adjusted R-squared: 0.003332
## F-statistic: 2.291 on 1 and 385 DF, p-value: 0.131
# Define the response variables
celltypes <- c("VSMC_I", "VSMC_II")
# Initialize a list to store the model results
macro_res_list <- list()
# Use a for loop to fit linear regression models for each response variable
for (ctypes in celltypes) {
formula <- paste(ctypes, '~ as.numeric([*])')
# Note: Replace 'as.numeric([*])' with the actual predictor variable
# you want to use. I've used a placeholder here.
macro_res <- RegParallel(
data = case2_aao,
formula = formula,
FUN = function(formula, data)
lm(formula = formula, data = data),
FUNtype = 'lm',
variables = colnames(case2_aao)[24:31],
blocksize = 8,
p.adjust = "BH",
cores = 2,
nestedParallel = FALSE,
conflevel = 95
)
# Store the results in the list
macro_res_list[[ctypes]] <- macro_res
}
##
## ##############################
## #RegParallel
## ##############################
## System is:
## -- Linux
## Blocksize:
## -- 8
## Cores / Threads:
## -- 2
## Terms included in model:
## -- VSMC_I
## First 5 formulae:
## -- VSMC_I ~ as.numeric(Pt_0.001)
## -- VSMC_I ~ as.numeric(Pt_0.05)
## -- VSMC_I ~ as.numeric(Pt_0.1)
## -- VSMC_I ~ as.numeric(Pt_0.2)
## -- VSMC_I ~ as.numeric(Pt_0.3)
## Processing 8 formulae, batch 1 of 1
## -- index1: 1; index2: 8
## Done!
##
## ##############################
## #RegParallel
## ##############################
## System is:
## -- Linux
## Blocksize:
## -- 8
## Cores / Threads:
## -- 2
## Terms included in model:
## -- VSMC_II
## First 5 formulae:
## -- VSMC_II ~ as.numeric(Pt_0.001)
## -- VSMC_II ~ as.numeric(Pt_0.05)
## -- VSMC_II ~ as.numeric(Pt_0.1)
## -- VSMC_II ~ as.numeric(Pt_0.2)
## -- VSMC_II ~ as.numeric(Pt_0.3)
## Processing 8 formulae, batch 1 of 1
## -- index1: 1; index2: 8
## Done!
PGS_aao_results <- dplyr::bind_rows(macro_res_list, .id ="celltypes")
PGS_aao_results
## celltypes Variable Term Beta StandardError
## <char> <char> <char> <num> <num>
## 1: VSMC_I Pt_0.001 as.numeric(Pt_0.001) -23.004867 29.02425
## 2: VSMC_I Pt_0.05 as.numeric(Pt_0.05) -329.177298 208.45690
## 3: VSMC_I Pt_0.1 as.numeric(Pt_0.1) -335.667745 305.66439
## 4: VSMC_I Pt_0.2 as.numeric(Pt_0.2) -494.007725 400.82453
## 5: VSMC_I Pt_0.3 as.numeric(Pt_0.3) -649.614535 476.74683
## 6: VSMC_I Pt_0.4 as.numeric(Pt_0.4) -737.089294 550.20681
## 7: VSMC_I Pt_0.5 as.numeric(Pt_0.5) -791.208475 611.74571
## 8: VSMC_I Pt_1 as.numeric(Pt_1) -1227.894285 811.31542
## 9: VSMC_II Pt_0.001 as.numeric(Pt_0.001) -5.092318 21.21451
## 10: VSMC_II Pt_0.05 as.numeric(Pt_0.05) -16.609081 152.74319
## 11: VSMC_II Pt_0.1 as.numeric(Pt_0.1) -108.040075 223.53348
## 12: VSMC_II Pt_0.2 as.numeric(Pt_0.2) -74.219052 293.30778
## 13: VSMC_II Pt_0.3 as.numeric(Pt_0.3) -88.330506 349.01711
## 14: VSMC_II Pt_0.4 as.numeric(Pt_0.4) 76.733855 402.77798
## 15: VSMC_II Pt_0.5 as.numeric(Pt_0.5) 133.929869 447.72593
## 16: VSMC_II Pt_1 as.numeric(Pt_1) 186.362007 594.25481
## t P OR ORlower ORupper P.adjust
## <num> <num> <num> <num> <num> <num>
## 1: -0.7926084 0.4284941 1.021206e-10 2.011987e-35 5.183242e+14 0.4284941
## 2: -1.5791144 0.1151310 1.096771e-143 3.992050e-321 3.012674e+34 0.2913632
## 3: -1.0981579 0.2728219 1.664756e-146 0.000000e+00 2.531319e+114 0.3117964
## 4: -1.2324788 0.2185224 2.852141e-215 0.000000e+00 NA 0.2913632
## 5: -1.3625986 0.1738056 7.516089e-283 0.000000e+00 NA 0.2913632
## 6: -1.3396586 0.1811466 7.692602e-321 0.000000e+00 NA 0.2913632
## 7: -1.2933617 0.1966615 0.000000e+00 0.000000e+00 NA 0.2913632
## 8: -1.5134610 0.1309830 0.000000e+00 0.000000e+00 NA 0.2913632
## 9: -0.2400394 0.8104274 6.143763e-03 5.377916e-21 7.018671e+15 0.9134665
## 10: -0.1087386 0.9134665 6.120232e-08 5.909207e-138 6.338793e+122 0.9134665
## 11: -0.4833284 0.6291373 1.198924e-47 6.408385e-238 2.243029e+143 0.9134665
## 12: -0.2530415 0.8003712 5.848915e-33 1.267767e-282 2.698430e+217 0.9134665
## 13: -0.2530836 0.8003387 4.350593e-39 0.000000e+00 5.277738e+258 0.9134665
## 14: 0.1905115 0.8490087 2.113925e+33 0.000000e+00 NA 0.9134665
## 15: 0.2991336 0.7649995 1.462188e+58 0.000000e+00 NA 0.9134665
## 16: 0.3136062 0.7539900 8.629615e+80 0.000000e+00 NA 0.9134665
##Now, the same above but with PGS of DAo distensibility
pgs_dao <- read.table("dis_pgs/dao.all_score", header = T)
cell_prop <- read.csv("/home/mchopra/Documents/PhD-Year1/deconvolution/Deconvolution_results/results/CIBERSORTx_Results.csv")
rownames(pgs_dao) <- pgs_dao[,1]
rownames(cell_prop) <- cell_prop[,1]
intersect <- intersect(rownames(cell_prop), rownames(pgs_dao))
cell_prop <- cell_prop[intersect, ]
pgs_dao <- pgs_dao[intersect,]
#changing the column1 name
names(pgs_dao)[1] <- "Mixture"
case2_dao <- merge(cell_prop, pgs_dao, by = "Mixture", all = TRUE)
#colnames(case2_dao)
# Define the response variables
celltypes <- c("VSMC_I", "VSMC_II")
# Initialize a list to store the model results
macro_res_list <- list()
# Use a for loop to fit linear regression models for each response variable
for (ctypes in celltypes) {
formula <- paste(ctypes, '~ as.numeric([*])')
# Note: Replace 'as.numeric([*])' with the actual predictor variable
# you want to use. I've used a placeholder here.
macro_res <- RegParallel(
data = case2_dao,
formula = formula,
FUN = function(formula, data)
lm(formula = formula, data = data),
FUNtype = 'lm',
variables = colnames(case2_aao)[24:31],
blocksize = 8,
p.adjust = "BH",
cores = 2,
nestedParallel = FALSE,
conflevel = 95
)
# Store the results in the list
macro_res_list[[ctypes]] <- macro_res
}
##
## ##############################
## #RegParallel
## ##############################
## System is:
## -- Linux
## Blocksize:
## -- 8
## Cores / Threads:
## -- 2
## Terms included in model:
## -- VSMC_I
## First 5 formulae:
## -- VSMC_I ~ as.numeric(Pt_0.001)
## -- VSMC_I ~ as.numeric(Pt_0.05)
## -- VSMC_I ~ as.numeric(Pt_0.1)
## -- VSMC_I ~ as.numeric(Pt_0.2)
## -- VSMC_I ~ as.numeric(Pt_0.3)
## Processing 8 formulae, batch 1 of 1
## -- index1: 1; index2: 8
## Done!
##
## ##############################
## #RegParallel
## ##############################
## System is:
## -- Linux
## Blocksize:
## -- 8
## Cores / Threads:
## -- 2
## Terms included in model:
## -- VSMC_II
## First 5 formulae:
## -- VSMC_II ~ as.numeric(Pt_0.001)
## -- VSMC_II ~ as.numeric(Pt_0.05)
## -- VSMC_II ~ as.numeric(Pt_0.1)
## -- VSMC_II ~ as.numeric(Pt_0.2)
## -- VSMC_II ~ as.numeric(Pt_0.3)
## Processing 8 formulae, batch 1 of 1
## -- index1: 1; index2: 8
## Done!
PGS_dao_results <- dplyr::bind_rows(macro_res_list, .id ="celltypes")
PGS_dao_results
## celltypes Variable Term Beta StandardError
## <char> <char> <char> <num> <num>
## 1: VSMC_I Pt_0.001 as.numeric(Pt_0.001) 47.185620 30.74294
## 2: VSMC_I Pt_0.05 as.numeric(Pt_0.05) 61.641115 228.71367
## 3: VSMC_I Pt_0.1 as.numeric(Pt_0.1) 122.438727 320.83086
## 4: VSMC_I Pt_0.2 as.numeric(Pt_0.2) 197.473880 456.60313
## 5: VSMC_I Pt_0.3 as.numeric(Pt_0.3) -136.908191 552.22005
## 6: VSMC_I Pt_0.4 as.numeric(Pt_0.4) -309.183116 631.54397
## 7: VSMC_I Pt_0.5 as.numeric(Pt_0.5) -359.923836 723.33861
## 8: VSMC_I Pt_1 as.numeric(Pt_1) -619.771979 987.69885
## 9: VSMC_II Pt_0.001 as.numeric(Pt_0.001) -31.213829 22.46645
## 10: VSMC_II Pt_0.05 as.numeric(Pt_0.05) -53.288117 167.04214
## 11: VSMC_II Pt_0.1 as.numeric(Pt_0.1) -18.493419 234.37166
## 12: VSMC_II Pt_0.2 as.numeric(Pt_0.2) -151.469693 333.48659
## 13: VSMC_II Pt_0.3 as.numeric(Pt_0.3) -179.762443 403.25998
## 14: VSMC_II Pt_0.4 as.numeric(Pt_0.4) -2.856055 461.41218
## 15: VSMC_II Pt_0.5 as.numeric(Pt_0.5) 73.441313 528.47044
## 16: VSMC_II Pt_1 as.numeric(Pt_1) 204.669520 721.69135
## t P OR ORlower ORupper P.adjust
## <num> <num> <num> <num> <num> <num>
## 1: 1.534843905 0.1256435 3.107811e+20 2.108692e-06 4.580322e+46 0.8043260
## 2: 0.269512160 0.7876800 5.893809e+26 1.227349e-168 2.830246e+221 0.8043260
## 3: 0.381630146 0.7029463 1.494389e+53 1.209880e-220 NA 0.8043260
## 4: 0.432484723 0.6656313 5.778516e+85 0.000000e+00 NA 0.8043260
## 5: -0.247923254 0.8043260 3.479589e-60 0.000000e+00 NA 0.8043260
## 6: -0.489567048 0.6247191 5.290280e-135 0.000000e+00 NA 0.8043260
## 7: -0.497586927 0.6190593 4.864789e-157 0.000000e+00 NA 0.8043260
## 8: -0.627490840 0.5307095 6.861983e-270 0.000000e+00 NA 0.8043260
## 9: -1.389353050 0.1655284 2.779754e-14 2.091850e-33 3.693874e+05 0.9950645
## 10: -0.319010027 0.7498919 7.198881e-24 4.685331e-166 1.106088e+119 0.9950645
## 11: -0.078906378 0.9371481 9.298446e-09 2.957286e-208 2.923664e+191 0.9950645
## 12: -0.454200249 0.6499407 1.650245e-66 0.000000e+00 1.207395e+218 0.9950645
## 13: -0.445773076 0.6560117 8.514570e-79 0.000000e+00 NA 0.9950645
## 14: -0.006189812 0.9950645 5.749514e-02 0.000000e+00 NA 0.9950645
## 15: 0.138969576 0.8895469 7.855194e+31 0.000000e+00 NA 0.9950645
## 16: 0.283597027 0.7768717 7.706254e+88 0.000000e+00 NA 0.9950645