How to create a (total) Ghrelin calibration curve?
- Version 13
- by (unknown)
- Version 14
- by (unknown)
Deletions or items before changed
Additions or items after changed
1 | == Computing Ghrelin calibration curve using R == | |||
---|---|---|---|---|
2 | ||||
3 | The calibration curve can be created using R, and libracy "drc". Optionally, one can use library "sfmisc" for formatting of the labels on plot axis. | |||
4 | + | |||
5 | + | We will assume a 4-parameter log-logistic model:
|
||
6 | + | [[Image(LL4.png)]]
|
||
7 | + | |||
8 | The R code is attached below. | |||
9 | ||||
10 | The example assumes the data to be available in file "ghrelin_conc std_a std_b avg.csv" | |||
11 | ||||
12 | The measured data: | |||
13 | ||||
14 | || '''Ghrelin (ng/ml)''' || '''Standard a''' || '''Standard b''' || | |||
15 | || 1000000 ||-0.040596823||-0.052699697 || | |||
16 | || 100000 ||0.136105144||0.119766263 || | |||
17 | || 10000 ||0.61356354||0.606906959 || | |||
18 | || 1000 ||0.846543873||0.839887292 || | |||
19 | || 100 ||0.887693646||0.88345764 || | |||
20 | ||0||0.896770802||0.896165658 || | |||
21 | ||||
22 | {{{ | |||
23 | ||||
24 | ##### Install libraries | |||
25 | install.packages("drc") | |||
26 | install.packages("sfsmisc") | |||
27 | require(drc) | |||
28 | library(sfsmisc) | |||
29 | ||||
30 | ##### Read the data | |||
31 | hormone.data <- read.csv("ghrelin_conc std_a std_b avg.csv") | |||
32 | hormone.data <- hormone.data[,1:3] | |||
33 | colnames(hormone.data)[1:3] <- c("Concentration","Response_1", "Response_2") | |||
34 | ||||
35 | ##### Reorganize the data | |||
36 | hormone.data <- reshape(hormone.data, varying=c("Response_1","Response_2"), direction="long", v.names=c("Response")) | |||
37 | hormone.data <- hormone.data[,c("Concentration", "Response")] | |||
38 | ||||
39 | ##### Fitting the model (4-parameter log-logistic function) | |||
40 | hormone.data.model <- drm(Response ~ Concentration, data = hormone.data, fct = LL.4()) | |||
41 | summary(hormone.data.model) | |||
42 | ||||
43 | }}} | |||
44 | ||||
45 | - | + | The resultant parameters of a log-logistic equation are:
|
|
46 | - | The resultant parameters of a log-logistic |
+ | |
47 | ||||
48 | {{{ | |||
49 | Model fitted: Log-logistic (ED50 as parameter) (4 parms) | |||
50 | Parameter estimates: | |||
51 | Estimate Std. Error t-value p-value | |||
52 | b:(Intercept) 9.5057e-01 2.2294e-02 4.2638e+01 0 | |||
53 | c:(Intercept) -7.6010e-02 6.9075e-03 -1.1004e+01 0 | |||
54 | d:(Intercept) 8.9163e-01 3.3216e-03 2.6843e+02 0 | |||
55 | e:(Intercept) 2.5221e+04 7.7727e+02 3.2448e+01 0 | |||
56 | }}} | |||
57 | ||||
58 | The calibration curve can be plotted using the commands below: | |||
59 | ||||
60 | {{{ | |||
61 | ##### Plotting a nice plot | |||
62 | par(pty="s", mar=c(5,5,1,1)) | |||
63 | plot(hormone.data.model, type="confidence", cex.lab=2, axes=F, xlim=c(-10,10^6)) | |||
64 | axis(side=1, at=hormone.data[1:6,1], labels=pretty10exp(hormone.data[1:6,1]), cex.axis=1.2) | |||
65 | axis(side=2, at=seq(0,1,0.2), labels=seq(0,1,0.2)) | |||
66 | plot(hormone.data.model, type="all", add=T, pch=21, col="red", lwd=1, cex=2, bg="green") | |||
67 | }}} | |||
68 | ||||
69 | [[Image(Ghrelin.png)]] | |||
70 | ||||
71 | - | The parameters of the eqution can be plugged into the formula below, and used in Excel, or other spreadsheet program.
|
+ | The parameters of the eqution can be plugged into the formula below (an inverse of the model), and used in Excel, or other spreadsheet program.
|
72 | - | + | [[Image(LL4-inv.png)]]
|
|
73 | ||||
74 | However, the concentration can be also easily estimated in R using "ED" function of the "drc" library. The code below demonstrates the concentration estimated from the response of 0.1, assuming alpha=0.05. The code returns the estimation, the error, and the condfidence interval. | |||
75 | ||||
76 | {{{ | |||
77 | ##### Computing the concentration from the response, for instance for a response=0.1, and alpha=1-0.95 | |||
78 | ED(hormone.data.model, respLev=0.1, interval="delta", type="absolute", level=0.95) | |||
79 | }}} |