2013-04-05 198 views
0

我已经写了一个简短的脚本来绘制2个独立的测量设备上测量的放射性活动的趋势。该脚本如下所示ggplot2 - 因子颜色和图例调整

pkgLoad <- function(x) 
    { 
    if (!require(x,character.only = TRUE)) 
    { 
     install.packages(x,dep=TRUE, repos='http://star-www.st-andrews.ac.uk/cran/') 
     if(!require(x,character.only = TRUE)) stop("Package not found") 
    } 

    } 

pkgLoad("ggplot2") 
pkgLoad("XLConnect") 
pkgLoad("reshape2") 

#Load the workbook 
wb<-loadWorkbook("CapintecQC.xlsx") 

df_blue <-readWorksheet(wb, sheet = "Blue", startCol=1, endCol=6) 
#sort date format 
df_blue$Date <- as.Date(df_blue$Date , "%d/%m/%y") 
df_blue[order(df_blue$Date),] 

df_gold <-readWorksheet(wb, sheet = "Gold", startCol=1, endCol=6) 
df_gold$Date <- as.Date(df_gold$Date , "%d/%m/%y") 
df_gold[order(df_gold$Date),] 

#Reference Cs-137 details 
ref_activity <- 9.3 
half_life <- 30.23 
ref_date <- as.Date('06/01/08',format='%d/%m/%y') 

blue_melt <- melt(df_blue[,c(1,2:6)], id="Date", value.name="Activity", variable.name="Isotope") 

#Add new column to data frame with expected activity 
df_gold["Exp_Act"] <- round(ref_activity*exp((-0.693/half_life)*as.numeric(difftime(df_gold$Date,ref_date))/365.25),3) 
df_gold["Exp_Act_0.95"] <- 0.95 * df_gold$Exp_Act 
df_gold["Exp_Act_1.05"] <- 1.05 * df_gold$Exp_Act 


gold_melt <- melt(df_gold[,c(1,2:6)], id="Date", value.name="Activity", variable.name="Isotope") 

p <- ggplot(NULL)+geom_point(data = gold_melt, aes(x=Date,y=Activity, col=Isotope)) + geom_ribbon(data = df_gold, aes(x = Date, ymin = Exp_Act_0.95, ymax = Exp_Act_1.05), fill='blue', alpha=0.2) + geom_point(data = blue_melt, aes(x=Date,y=Activity, col=Isotope), shape=2) + theme_bw() 
print(p) 

enter image description here

我不是很称职的R/GGPLOT2。我希望最后一幅图显示每种放射性核素的测量活性对于这两种装置是相同的颜色(即红色的Cs-137,蓝色的99mTc)。我怎样才能做到这一点,因为我的图形绘制了不同的颜色。

此外,传说是不愉快的。 (i)从excel表头获取的每个核素的格式从Cs-137变为Cs.137。我怎样才能有Cs-137,Tc-99m等标题? (ii)每个放射性核素都在图例中复制 - 每个设备一个。是否有可能只显示第一个数据帧(df_gold)的传说或更好的只是文字的传说,在情节相匹配的标记颜色的文本颜色?)

df_gold结构

structure(list(Date = structure(c(15708, 15709, 15712, 15713, 
15714, 15715, 15716, 15719, 15720, 15721, 15722, 15723, 15726, 
15727, 15729, 15730, 15733, 15734, 15735, 15736, 15740, 15741, 
15743, 15747, 15748, 15749, 15750, 15751, 15754, 15755, 15756, 
15757, 15758, 15761, 15762, 15764, 15765, 15768, 15769, 15770, 
15771, 15772, 15775, 15776, 15777, 15779, 15782, 15783, 15784, 
15785, 15786, 15789, 15790, 15791, 15792, 15797, 15798, 15799, 
15800), class = "Date"), Cs..137 = c(8.2, 8.1, 8.1, 8.1, 8.1, 
8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 
8.1, 8.1, 8.1, 8.2, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8, 8.2, 
8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 
8.1, 8.1, 8.1, 8, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 8.1, 
8.1, 8.1), In..111 = c(6.49, 6.47, 6.48, 6.43, 6.49, 6.51, 6.5, 
6.47, 6.48, 6.4, 6.48, 6.48, 6.48, 6.49, 6.49, 6.47, 6.48, 6.48, 
6.5, 6.47, 6.49, 6.55, 6.46, 6.49, 6.48, 6.48, 6.46, 6.48, 6.49, 
6.44, 6.49, 6.46, 6.45, 6.46, 6.46, 6.43, 6.49, 6.47, 6.45, 6.43, 
6.44, 6.44, 6.44, 6.46, 6.45, 6.47, 6.45, 6.43, 6.44, 6.47, 6.45, 
6.46, 6.45, 6.46, 6.39, 6.46, 6.44, 6.42, 6.41), I..123 = c(6.97, 
6.94, 6.96, 6.91, 6.92, 6.95, 6.93, 6.92, 6.93, 7, 6.97, 6.96, 
6.96, 6.94, 6.98, 6.97, 6.95, 6.95, 6.94, 6.96, 6.97, 7.01, 6.92, 
7, 6.98, 6.97, 6.91, 6.99, 6.95, 6.88, 6.96, 6.91, 6.91, 6.93, 
6.94, 6.94, 6.97, 6.93, 6.93, 6.93, 6.96, 6.94, 6.94, 6.92, 6.93, 
6.91, 6.93, 6.92, 6.92, 6.91, 6.91, 6.89, 6.92, 6.9, 6.9, 6.91, 
6.91, 6.9, 6.9), I..131 = c(10.5, 10.5, 10.5, 10.5, 10.5, 10.5, 
10.5, 10.5, 10.5, 10.8, 10.5, 10.6, 10.5, 10.5, 10.5, 10.5, 10.5, 
10.5, 10.5, 10.5, 10.5, 10.6, 10.5, 10.5, 10.5, 10.5, 10.5, 10.5, 
10.5, 10.5, 10.5, 10.5, 10.5, 10.5, 10.5, 10.5, 10.5, 10.5, 10.4, 
10.5, 10.4, 10.5, 10.5, 10.5, 10.4, 10.5, 10.4, 10.4, 10.5, 10.4, 
10.4, 10.4, 10.4, 10.4, 10.3, 10.5, 10.5, 10.5, 10.6), Tc..99m = c(15, 
15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15.1, 15, 15, 15.1, 15, 
15, 15, 15, 15.1, 15, 15.1, 15, 15, 15, 15, 15, 15, 15, 15, 15, 
15, 15, 15, 15, 14.9, 15, 15, 15, 15, 15, 15, 15, 15, 15, 15, 
14.9, 14.8, 14.9, 14.9, 14.9, 14.9, 15, 15, 14.8, 15, 15, 15, 
15), Exp_Act = c(8.294, 8.293, 8.292, 8.291, 8.291, 8.29, 8.29, 
8.288, 8.288, 8.287, 8.287, 8.286, 8.285, 8.284, 8.283, 8.283, 
8.281, 8.28, 8.28, 8.279, 8.277, 8.277, 8.276, 8.274, 8.273, 
8.273, 8.272, 8.272, 8.27, 8.27, 8.269, 8.269, 8.268, 8.266, 
8.266, 8.265, 8.264, 8.263, 8.262, 8.262, 8.261, 8.261, 8.259, 
8.259, 8.258, 8.257, 8.256, 8.255, 8.255, 8.254, 8.254, 8.252, 
8.251, 8.251, 8.25, 8.248, 8.247, 8.247, 8.246), Exp_Act_0.95 = c(7.8793, 
7.87835, 7.8774, 7.87645, 7.87645, 7.8755, 7.8755, 7.8736, 7.8736, 
7.87265, 7.87265, 7.8717, 7.87075, 7.8698, 7.86885, 7.86885, 
7.86695, 7.866, 7.866, 7.86505, 7.86315, 7.86315, 7.8622, 7.8603, 
7.85935, 7.85935, 7.8584, 7.8584, 7.8565, 7.8565, 7.85555, 7.85555, 
7.8546, 7.8527, 7.8527, 7.85175, 7.8508, 7.84985, 7.8489, 7.8489, 
7.84795, 7.84795, 7.84605, 7.84605, 7.8451, 7.84415, 7.8432, 
7.84225, 7.84225, 7.8413, 7.8413, 7.8394, 7.83845, 7.83845, 7.8375, 
7.8356, 7.83465, 7.83465, 7.8337), Exp_Act_1.05 = c(8.7087, 8.70765, 
8.7066, 8.70555, 8.70555, 8.7045, 8.7045, 8.7024, 8.7024, 8.70135, 
8.70135, 8.7003, 8.69925, 8.6982, 8.69715, 8.69715, 8.69505, 
8.694, 8.694, 8.69295, 8.69085, 8.69085, 8.6898, 8.6877, 8.68665, 
8.68665, 8.6856, 8.6856, 8.6835, 8.6835, 8.68245, 8.68245, 8.6814, 
8.6793, 8.6793, 8.67825, 8.6772, 8.67615, 8.6751, 8.6751, 8.67405, 
8.67405, 8.67195, 8.67195, 8.6709, 8.66985, 8.6688, 8.66775, 
8.66775, 8.6667, 8.6667, 8.6646, 8.66355, 8.66355, 8.6625, 8.6604, 
8.65935, 8.65935, 8.6583)), row.names = c(NA, -59L), .Names = c("Date", 
"Cs..137", "In..111", "I..123", "I..131", "Tc..99m", "Exp_Act", 
"Exp_Act_0.95", "Exp_Act_1.05"), class = "data.frame") 

df_blue结构

structure(list(Date = structure(c(15790, 15791, 15792, 15797, 
15798, 15799, 15800), class = "Date"), Cs.137 = c(8.1, 8.2, 8.2, 
8.2, 8.2, 8.2, 8.2), I.123 = c(6.82, 6.85, 6.91, 6.84, 6.82, 
6.82, 6.83), I.131 = c(10.5, 10.6, 10.6, 10.5, 10.6, 10.6, 10.6 
), In.111 = c(6.35, 6.45, 6.43, 6.37, 6.38, 6.4, 6.37), X99m.Tc = c(15, 
15, 15.1, 15.1, 15.1, 15.1, 15.1)), .Names = c("Date", "Cs.137", 
"I.123", "I.131", "In.111", "X99m.Tc"), row.names = c(NA, -7L 
), class = "data.frame") 

回答

1

我的方法对双方都是数据帧结合在一起,然后添加包含的装置(金或蓝色)名称新列。

df<-rbind(gold_melt,blue_melt) 
df$device<-rep(c("gold","blue"),c(nrow(gold_melt),nrow(blue_melt))) 

随着从Isotopecar更改名称功能recode(),因为他们应该。

df$Isotope<-recode(df$Isotope,"c('Cs..137','Cs.137')='Cs-137'; 
     c('I..123','I.123')='I-123'; 
     c('I..131','I.131')='I-131'; 
     c('In..111','In.111')='In-111' 
     ;c('Tc..99m','X99m.Tc')='Tc-99m'") 

现在你只需要使用新的数据帧调用geom_point()。我还添加了shape=device以获得每个设备的不同形状。

ggplot(NULL) + 
    geom_point(data=df,aes(x=Date,y=Activity, col=Isotope,shape=device))+ 
    geom_ribbon(data = df_gold, aes(x = Date, ymin = Exp_Act_0.95, ymax = Exp_Act_1.05), fill='blue', alpha=0.2) 

enter image description here

+0

感谢。有没有一个简单的方法来设置设备形状为黄金开放圈。你的回复很有启发性。 – moadeep 2013-04-05 12:11:07

+1

您可以通过添加+ scale_shape_manual(values = c(16,1))来更改形状 – 2013-04-05 12:14:30

+0

非常感谢 – moadeep 2013-04-05 12:16:42

1

万一你想 “导火索” 的传说一起,然后在Didzis的回答建立:

df <- transform(df, device = factor(device, levels=unique(device)), 
          grp = paste(Isotope, device, sep="_")) 

require(RColorBrewer) 
ggplot() + geom_point(data = df, aes(x = Date, y = Activity, 
         colour=grp, shape = grp, fill=grp)) + 
geom_ribbon(data = df_gold, aes(x = Date, ymin = Exp_Act_0.95, 
        ymax = Exp_Act_1.05), fill='blue', alpha=0.2) + 
scale_shape_manual("", values=rep(c(21,24), 5)) + 
scale_fill_manual("", values=rep(brewer.pal(5, "Set1"), each=2)) + 
scale_colour_manual("", values=rep(brewer.pal(5, "Set1"), each=2)) 

enter image description here