-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathregression.r
140 lines (112 loc) · 5.31 KB
/
regression.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
# """
# Perform multi-linear regression for all the run down rates with temperature,
# thold, and INaCa as the indepedent variables
# """
# get path
pathwork = getwd()
pathrel = "output/r_rate_database.csv"
path = paste(pathwork, pathrel, sep = "/")
# load data
dataset = read.csv(path)
dataset['thold'] = log(dataset['thold'])
# Reorder the levels to set 310 as the reference level
dataset$Temperature = factor(dataset$Temperature, levels = c(310, 298))
# perform multi-linear regression
multi.fit = lm(Run.rate~thold + factor(Temperature) + factor(INaCa), data=dataset)
res = summary(multi.fit)
# print output
print("############### Summary of model fit ###############")
print(summary(multi.fit))
print(confint(multi.fit))
print("####################################################")
# For H1: beta < 0
p_value_less = pt(coef(res)[, 3], multi.fit$df, lower = TRUE)
print("## Test for directionality (beta coefficients <0) ##")
print(p_value_less)
print("####################################################")
# Statistical reduction by addition of temperature
initial.model <- lm(Run.rate ~ 1, data = dataset)
full.model <- lm(Run.rate ~ factor(Temperature), data = dataset)
# Perform F-test
f_test_result <- anova(initial.model, full.model)
print("##### Statistical reduction due to temperature #####")
print(f_test_result)
print("####################################################")
# Statistical reduction by addition of INaCa
initial.model <- lm(Run.rate ~ factor(Temperature), data = dataset)
full.model <- lm(Run.rate ~ factor(Temperature) + factor(INaCa), data = dataset)
# Perform F-test
f_test_result <- anova(initial.model, full.model)
print("######## Statistical reduction due to INaCa ########")
print(f_test_result)
print("####################################################")
# Statistical reduction by addition of thold
initial.model <- lm(Run.rate ~ factor(Temperature)+ factor(INaCa), data = dataset)
full.model <- lm(Run.rate ~ factor(Temperature) + factor(INaCa) + thold, data = dataset)
# Perform F-test
f_test_result <- anova(initial.model, full.model)
print("######## Statistical reduction due to thold ########")
print(f_test_result)
print("####################################################")
# For H1: beta > 0
# p_value_greater = pt(coef(res)[, 3], multi.fit$df, lower = FALSE)
# print(p_value_greater)
#%%%%%%%%%%%%%%%%%%%%%%%%OUTPUT%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
# [1] "############### Summary of model fit ###############"
# Call:
# lm(formula = Run.rate ~ thold + factor(Temperature) + factor(INaCa),
# data = dataset)
# Residuals:
# Min 1Q Median 3Q Max
# -0.158729 -0.016681 0.003505 0.026190 0.121515
# Coefficients:
# Estimate Std. Error t value Pr(>|t|)
# (Intercept) 0.158795 0.019715 8.055 2.16e-13 ***
# thold -0.016373 0.006545 -2.502 0.01342 *
# factor(Temperature)298 -0.062213 0.007123 -8.734 4.16e-15 ***
# factor(INaCa)On -0.021030 0.007023 -2.994 0.00321 **
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# Residual standard error: 0.04355 on 152 degrees of freedom
# Multiple R-squared: 0.4045, Adjusted R-squared: 0.3928
# F-statistic: 34.42 on 3 and 152 DF, p-value: < 2.2e-16
# 2.5 % 97.5 %
# (Intercept) 0.11984411 0.197745292
# thold -0.02930387 -0.003442791
# factor(Temperature)298 -0.07628548 -0.048139628
# factor(INaCa)On -0.03490557 -0.007153578
# [1] "####################################################"
# [1] "## Test for directionality (beta coefficients <0) ##"
# (Intercept) thold factor(Temperature)298 factor(INaCa)On
# 1.000000e+00 6.708997e-03 2.079688e-15 1.606375e-03
# [1] "####################################################"
# [1] "##### Statistical reduction due to temperature #####"
# Analysis of Variance Table
# Model 1: Run.rate ~ 1
# Model 2: Run.rate ~ factor(Temperature)
# Res.Df RSS Df Sum of Sq F Pr(>F)
# 1 155 0.48405
# 2 154 0.31631 1 0.16773 81.661 6.389e-16 ***
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# [1] "####################################################"
# [1] "######## Statistical reduction due to INaCa ########"
# Analysis of Variance Table
# Model 1: Run.rate ~ factor(Temperature)
# Model 2: Run.rate ~ factor(Temperature) + factor(INaCa)
# Res.Df RSS Df Sum of Sq F Pr(>F)
# 1 154 0.31631
# 2 153 0.30010 1 0.016211 8.2647 0.004619 **
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# [1] "####################################################"
# [1] "######## Statistical reduction due to thold ########"
# Analysis of Variance Table
# Model 1: Run.rate ~ factor(Temperature) + factor(INaCa)
# Model 2: Run.rate ~ factor(Temperature) + factor(INaCa) + thold
# Res.Df RSS Df Sum of Sq F Pr(>F)
# 1 153 0.30010
# 2 152 0.28824 1 0.011868 6.2586 0.01342 *
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# [1] "####################################################"