How to create a (total) Ghrelin calibration curve?
 Version 21
 by (unknown)
 Version 22
 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    <math>n = 5</math>


6  
7  We will assume a 4parameter loglogistic model:  
8  
9  [[Image(LL4.png, 300px)]]  
10  
11  The R code is attached below.  
12  
13  The example assumes the data to be available in file "ghrelin_conc std_a std_b avg.csv"  
14  
15  The measured data:  
16  
17   '''Ghrelin (ng/ml)'''  '''Standard a'''  '''Standard b'''   
18   1000000 0.0405968230.052699697   
19   100000 0.1361051440.119766263   
20   10000 0.613563540.606906959   
21   1000 0.8465438730.839887292   
22   100 0.8876936460.88345764   
23  00.8967708020.896165658   
24  
25  {{{  
26  
27  ##### Install libraries  
28  install.packages("drc")  
29  install.packages("sfsmisc")  
30  require(drc)  
31  library(sfsmisc)  
32  
33  ##### Read the data  
34  hormone.data < read.csv("ghrelin_conc std_a std_b avg.csv")  
35  hormone.data < hormone.data[,1:3]  
36  colnames(hormone.data)[1:3] < c("Concentration","Response_1", "Response_2")  
37  
38  ##### Reorganize the data  
39  hormone.data < reshape(hormone.data, varying=c("Response_1","Response_2"), direction="long", v.names=c("Response"))  
40  hormone.data < hormone.data[,c("Concentration", "Response")]  
41  
42  ##### Fitting the model (4parameter loglogistic function)  
43  hormone.data.model < drm(Response ~ Concentration, data = hormone.data, fct = LL.4())  
44  summary(hormone.data.model)  
45  
46  }}}  
47  
48  The resultant parameters of a loglogistic equation are:  
49  
50  {{{  
51  Model fitted: Loglogistic (ED50 as parameter) (4 parms)  
52  Parameter estimates:  
53  Estimate Std. Error tvalue pvalue  
54  b:(Intercept) 9.5057e01 2.2294e02 4.2638e+01 0  
55  c:(Intercept) 7.6010e02 6.9075e03 1.1004e+01 0  
56  d:(Intercept) 8.9163e01 3.3216e03 2.6843e+02 0  
57  e:(Intercept) 2.5221e+04 7.7727e+02 3.2448e+01 0  
58  }}}  
59  
60  The calibration curve can be plotted using the commands below:  
61  
62  {{{  
63  ##### Plotting a nice plot  
64  par(pty="s", mar=c(5,5,1,1))  
65  plot(hormone.data.model, type="confidence", cex.lab=2, axes=F, xlim=c(10,10^6))  
66  axis(side=1, at=hormone.data[1:6,1], labels=pretty10exp(hormone.data[1:6,1]), cex.axis=1.2)  
67  axis(side=2, at=seq(0,1,0.2), labels=seq(0,1,0.2))  
68  plot(hormone.data.model, type="all", add=T, pch=21, col="red", lwd=1, cex=2, bg="green")  
69  }}}  
70  
71  [[Image(Ghrelin.png)]]  
72  
73  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.  
74  
75  [[Image(LL4inv.png, 230px)]]  
76  
77  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.  
78  
79  {{{  
80  ##### Computing the concentration from the response, for instance for a response=0.1, and alpha=10.95  
81  ED(hormone.data.model, respLev=0.1, interval="delta", type="absolute", level=0.95)  
82  }}} 