s-news
[Top] [All Lists]
<prev [Date] next>
[Advanced]
 <prev [Thread] next>

Re: trellis line plots with error bars

from [Jordan, Paul] [Permanent Link][Original]
To: "Thompson, David (MNR)" <david.thompson@mnr.gov.on.ca>, "s-news" <s-news@lists.biostat.wustl.edu>
Subject: Re: trellis line plots with error bars
From: "Jordan, Paul" <paul.jordan@roche.com>
Date: Thu, 9 Mar 2006 20:58:37 +0100
Thread-index: AcY9bInsMAWuHVt8S3Gm25MLQnS2EAGRoqoA
Thread-topic: trellis line plots with error bars 
Hi,

I propose to use an own-written panel function, such as:

"panel.MeanErrorBar"<-
function(x, y, type = "SD", level = 0.95, colLine = 1, ...)
{
        x1 <- ordered(x)
        if(type == "SD") {
                my <- tapply(y, x1, mean, na.rm = T)
                y.low <- y.upp <- tapply(y, x1, stdev, na.rm = T)
                Iy.low <- my - y.low
                Iy.upp <- my + y.low
        }
        if(type == "SEM") {
                my <- tapply(y, x1, mean, na.rm = T)
                N <- tapply(y, x1, function(z)
                length(z[!is.na(z)]))
                print(my)
                print(N)
                y.low <- y.upp <- tapply(y, x1, stdev, na.rm =
T)/sqrt(N)
                Iy.low <- my - y.low
                Iy.upp <- my + y.low
        }
        if(type == "CI") {
                alpha <- 1 - (1 - level)/2
                my <- tapply(y, x1, mean, na.rm = T)
                N <- tapply(y, x1, function(z)
                length(z[!is.na(z)]))
                y.low <- y.upp <- (qt(alpha, N - 1) * tapply(y, x1,
stdev,
                        na.rm = T))/sqrt(N)
                Iy.low <- my - y.low
                Iy.upp <- my + y.low
        }
        if(type == "CIlogn") {
                alpha <- 1 - (1 - level)/2
                y1 <- y[!is.na(y) & y > 0]
                xmy <- tapply(log(y1), x1, mean, na.rm = T)
                N <- tapply(y, x1, function(z)
                length(z[!is.na(z)]))
                y.low <- y.upp <- (qt(alpha, N - 1) * tapply(log(y1),
x1, stdev,
                        na.rm = T))/sqrt(N)
                my <- exp(xmy)
                Iy.low <- exp(xmy - y.low)
                Iy.upp <- exp(xmy + y.low)
        }
        xx <- unique(x[!is.na(x)])
        points(xx, my, ...)
        for(i in 1:length(xx)) {
                segments(xx[i], Iy.low[i], xx[i], Iy.upp[i], col =
colLine)
        }
}


An example of the use:

xyplot(AUC~Dose|diet
                ,data=PKdata
                ,panel=function(x,y)
                {
                        panel.xyplot(x,y,col=5)
                        panel.MeanErrorBar(x,y,type="CIlogn")
                }
                ,main="Test for panel.MeanErrorBar"
        )
The corresponding help file is attached.

Have fun!

Paul


-----Original Message-----
From: Thompson, David (MNR) [mailto:david.thompson@mnr.gov.on.ca] 
Sent: Wednesday, March 01, 2006 9:13 PM
To: s-news
Subject: trellis line plots with error bars


Hello,

I am wondering how to convert a tedious plotting routine to a trellis
method. The data has six columns of treatment means, six columns of
standard errors, one integer column 'year' (4 levels) and one factor
column 'group' (4 levels). (Data 'dump' below sig.)

At this time I, unfortunately, work mostly through the GUI and would use
a routine such as:

1) with a single species group subset,
   individually select (control+click) year, then mean1, then stderr1,
2) select 'Error Bar' plot from the '2D Plots' toolbar,
3) select the graph generated,
4) individually select year, then mean2, then stderr2,
5) shift+click on 'Error Bar' plot button.
6) repeat 4-5 until all six mean-stderr pairs have been plotted.

7) repeat 1-6 until all four species group subsets have been processed.

8) arrange the four plots on a single page.

'AARRGH' I *know* there have to be better ways!

Any and all help will be *much* appreciated.

Thanks, DaveT.
**********************************************************
Silviculture Data Analyst
Ontario Forest Research Institute
Ontario Ministry of Natural Resources
Sault Ste. Marie, Ontario, Canada
david.thompson@mnr.gov.on.ca
**********************************************************

P.S. Why does the output of data.dump go to my original user directory
rather than my current working directory (chapter)?

The contents of data.dump:
## Dump S Version 4 Dump ##
test
.named_I
structure
2
data.frame
structure
5
.named_I
structure
2
list
list
14
numeric
numeric
16
54.666666666666664
55.444444444444443
57.444444444444443
57.888888888888886
24
24.444444444444443
25.666666666666668
26
16.222222222222221
16.333333333333332
16.444444444444443
16.444444444444443
14.444444444444445
14.666666666666666
15.333333333333334
15.444444444444445
numeric
numeric
16
40
40.200000000000003
44.600000000000001
45.200000000000003
14.6
14.6
17
17
12.800000000000001
13
13.6
13.6
12.6
12.6
14
14.6
numeric
numeric
16
39.399999999999999
39.799999999999997
37
37.399999999999999
15.800000000000001
15.800000000000001
14.199999999999999
14.6
12.6
13.199999999999999
11.4
11.4
11
10.800000000000001
11.4
11.4
numeric
numeric
16
34.111111111111114
34.666666666666664
37.777777777777779
38.333333333333336
11.111111111111111
11
12.555555555555555
13
12.666666666666666
13.333333333333334
14.222222222222221
14.222222222222221
10.333333333333334
10.333333333333334
11
11.111111111111111
numeric
numeric
16
22.222222222222221
22.111111111111111
25.666666666666668
26.888888888888889
7
6.8888888888888893
8.4444444444444446
8.8888888888888893
7.8888888888888893
7.8888888888888893
8.8888888888888893
9.5555555555555554
7.333333333333333
7.333333333333333
8.3333333333333339
8.4444444444444446
numeric
numeric
16
30
30.222222222222221
32.444444444444443
33.333333333333336
12
12
12.888888888888889
13.222222222222221
9
9.1111111111111107
9.4444444444444446
9.5555555555555554
9
9.1111111111111107
10.111111111111111
10.555555555555555
integer
integer
16
1996
1997
1998
2000
1996
1997
1998
2000
1996
1997
1998
2000
1996
1997
1998
2000
numeric
numeric
16
2.1343747458109497
2.0419973277210119
2.2919612605262523
2.2939803135625088
1.6749792701868149
1.625415426480866
1.4433756729740645
1.4907119849998598
1.4021588293471072
1.4142135623730949
1.405456737852613
1.3955396318978297
0.83518313213182827
0.74535599249992979
0.8498365855987976
0.86780551954518381
numeric
numeric
16
2.9495762407505253
3.1208973068654471
3.0099833886584824
3.2310988842807022
2.1587033144922905
2.1587033144922905
1.7606816861659009
1.7606816861659009
0.37416573867739411
0.54772255750516607
1.1224972160321824
1.1224972160321824
1.0295630140987
1.0295630140987
0.63245553203367588
0.97979589711327109
numeric
numeric
16
2.9427877939124323
3.2771939216347876
1.51657508881031
1.16619037896906
1.6552945357246849
1.6552945357246849
1.3190905958272918
1.2083045973594571
0.92736184954957024
1.1575836902790224
0.74833147735478822
0.74833147735478822
1.2247448713915889
1.1575836902790224
0.6782329983125267
0.92736184954957024
numeric
numeric
16
1.3170299545699664
1.3437096247164249
1.4508405396889552
1.6666666666666667
1.2521586299538494
1.2360330811826106
0.91455986850857374
0.89752746785575066
0.57735026918962562
0.57735026918962573
0.61864048475889133
0.66202084932294369
1
1
1.0408329997330663
1.0467468711898207
numeric
numeric
16
0.982878111830785
0.99225394997176408
1.2247448713915892
0.82402205412174023
0.81649658092772592
0.82402205412174057
0.66897747659957207
0.65499034014175528
0.48432210483785254
0.48432210483785254
0.78959280019732747
0.60348780506667865
0.57735026918962584
0.57735026918962584
0.23570226039551587
0.29397236789606557
numeric
numeric
16
2.6977356760397746
2.7017369812861838
2.2798093920759102
2.2360679774997898
1.7400510848184252
1.7400510848184252
1.5041095967423201
1.4886401118860293
0.78173595997057166
0.77180244385832253
0.72860428047800008
0.76578048622723449
0.81649658092772592
0.80699105812966387
0.58794473579213113
0.50307695211874537
.named_I
structure
2
factor
structure
3
integer
integer
16
1
1
1
1
2
2
2
2
3
3
3
3
4
4
4
4
character
character
4
all
cryp
herb
wood
character
character
1
factor
character
character
3
.Data
.Label
class
character
character
14
value
value
value
value
value
value
value
value
value
value
value
value
value
value
character
character
14
CN
WR
SM
BM
SC2
SC1
yr
CNse
WRse
SMse
BMse
SC2se
SC1se
grp
character
character
16
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
list
list
14
character
character
5
TXPROP_ColWidth<-7
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-7
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-7
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-7
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-7
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-7
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
3
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
character
character
5
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
5
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Right
TXPROP_ObjectNote<-
TXPROP_ColFloatFormat<-Decimal
TXPROP_ColPrecision<-2
character
character
3
TXPROP_ColWidth<-6
TXPROP_ColJustification<-Left
TXPROP_ObjectNote<-
character
character
1
data.frame
character
character
5
.Data
.Names
row.names
.guiColInfo
class
panel.MeanErrorBar
Language Reference for PHM functions

Adding Means and Error Bars to a Panel

DESCRIPTION:

USAGE: 

panel.MeanErrorBar(x, y, type="SD", level=0.95, colLine=1, ...)

REQUIRED ARGUMENTS:

x
horizontal coordinates of data points.
y
vertical coordinates of data points.

OPTIONAL ARGUMENTS:

type
Type of the error bars. Currently the following types are available:
  1. "SD": Mean +- SD
  2. "SEM": Mean +- Standard Error of the Mean
  3. "CI": Mean and level*100% Confidence Interval for the Normal Distribution
  4. "CIlogn": Geometric Mean and level*100% Confidence Interval for the Lognormal Distribution
level
an optional numeric value for the interval confidence level. Defaults to 0.95.
colLine
colour number (see par for details). Defaults to black.
...
graphical parameters that are passed to points or lines.

SIDE EFFECTS:

The panel.MeanErrorBar function draws symbols for means and vertictal error bars on the display.

EXAMPLES: 

PKdata             <- PK
PKdata$Dose     <- as.numeric(as.character(PK$dose))
PKdata$AUC      <- exp(PK$log.AUC)

xyplot(AUC~Dose|diet
                ,data=""
                ,panel=function(x,y)
                {
                        panel.xyplot(x,y,col=5)
                        panel.MeanErrorBar(x,y,type="CIlogn")
                }
                ,main="Test for panel.MeanErrorBar"
        )
[More with this subject...]
<Prev in Thread] Current Thread [Next in Thread>
Previous by Date:  Biomedical Track, Schedule and Highlights, San Diego Data Mining Conference, March 29-31, San Diego, California, Lisa Solomon
Next by Date:  SPlus 7 Workbench and Eclipse 3.1, stephen . veronneau
Previous by Thread:  trellis line plots with error bars, Thompson, David (MNR)
Next by Thread:  SPlus 7 Workbench and Eclipse 3.1, stephen . veronneau
Indexes:  [Date] [Thread] [Top] [All Lists]