Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
380 views
in Technique[技术] by (71.8m points)

R script : add a padding to the ymax of the plot with ggplot2

I am trying to add a padding of 20 for the maximum of 5 curves I am plotting.

library(ggplot2)

... some code to compute 5 (bins) distributions y3(nSample,bin)

color_curve <- c("red", "green", "grey", "black")

max_y <- as.double(which.max(density(y3))+20)

# Save PDF - Chisquare
pdf(file="Chisquare_Distribution.pdf")
for (i in 1:nRed) {
  if (i == 1)
    plot(density(y3[,i]), col="blue", xlim=c(min(y3),max(y3)), ylim=c(0,max_y), main="z= 0.9595, 1.087, 1.2395, 1.45, 1.688")
  else
    lines(density(y3[,i]), col=color_curve[i-1])
}
dev.off()

But I get the following error at execution :

Error in which.max(density(y3)) :
  'list' object cannot be coerced to type 'double'
Execution halted

I would like to add this padding of 20 to the maximum of all 5 distributions but it fails with this error.

How to fix this ?

Update 1

Thanks for the suggested answer. Unfortunately, now, I get the following figure :

overestimated y_max limit of the plot

As you can see, the maximum in y_limit is too high for the 5 distributions, I don't know where it could come from.

Update 2

With the new command, I have the following figure :

under estimated of y_max of distribution

I get an under-estimated value for searching the max among the 5 distributions.

Edit

I provide the entire code to generate the plots with the input file:

my_data <- read.delim("Array_total_WITH_Shot_Noise.txt", header = FALSE, sep = " ")
array_2D <- array(my_data)
z_ph <- c(0.9595, 1.087, 1.2395, 1.45, 1.688)
b_sp <- c(1.42904922, 1.52601862, 1.63866958, 1.78259615, 1.91956918)
b_ph <- c(sqrt(1+z_ph))
ratio_squared <- (b_sp/b_ph)^2

nRed <- 5
nRow <- NROW(my_data)

nSample_var <- 1000000
nSample_mc <- 1000

Cl<-my_data[,2:length(my_data)]#suppose cl=var(alm)
Cl_sp <- array(0, dim=c(nRow,nRed))
Cl_ph <- array(0, dim=c(nRow,nRed))
length(Cl)
for (i in 1:length(Cl)) {
  #(shape/rate) convention : 
  Cl_sp[,i] <-(Cl[, i] * ratio_squared[i])
  Cl_ph[,i] <- (Cl[, i])
}
L <- array_2D[,1]
L <- 2*(array_2D[,1])+1

# Weighted sum of Chi squared distribution
y3_1<-array(0,dim=c(nSample_var,nRed));y3_2<-array(0,dim=c(nSample_var,nRed));y3<-array(0,dim=c(nSample_var,nRed));
  for (i in 1:nRed) {
    for (j in 1:nRow) {
      # Try to summing all the random variable
      y3_1[,i] <- y3_1[,i] + Cl_sp[j,i] * rchisq(nSample_var,df=L[j])
      y3_2[,i] <- y3_2[,i] + Cl_ph[j,i] * rchisq(nSample_var,df=L[j])
    }
    y3[,i] <- y3_1[,i]/y3_2[,i]
  }

print(paste0('n=',nSample_mc*nSample_var))
for (i in 1:nRed) {
  # compute the standard deviation of the ratio by Monte-Carlo
  print(paste0('red=',i,',mean_fid = ', ratio_squared[i],',mean_exp = ', mean(y3[,i])))
  print(paste0('numerator : var = ', var(y3_1[,i]), ', sigma = ', sd(y3_1[,i])))
  print(paste0('denominator : var = ', var(y3_2[,i]), ', sigma = ', sd(y3_2[,i])))
  print(paste0('var = ', var(y3[,i]), ', sigma = ', sd(y3[,i])))
}
print('#############################################################')
# par(mfrow=c(2,nRed))
color_curve <- c("red", "green", "grey", "black")

# Save PDF - Chisquare
pdf(file="Chisquare_Distribution.pdf")
for (i in 1:nRed) {
  if (i == 1) 
    plot(density(y3[,i]), col="blue", xlim=c(min(y3),max(y3)), main="z= 0.9595, 1.087, 1.2395, 1.45, 1.688")
  else 
    lines(density(y3[,i]), col=color_curve[i-1])
}
dev.off()

Hoping this will help you to do a version ggplot2 of the classical R plots.

See Question&Answers more detail:os

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Reply

0 votes
by (71.8m points)

No code to generate your sample data was provided, so I used the code @akrun had used previously: y3 <- matrix(rnorm(5000), ncol = 5)

library(tidyverse)
as.data.frame(y3) %>%
    mutate(row = row_number()) %>%     # add row to simplify next step
    pivot_longer(-row) %>%             # reshape long
    ggplot(aes(value, color = name)) + # map x to value, color to name
    geom_density()

enter image description here


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
OGeek|极客中国-欢迎来到极客的世界,一个免费开放的程序员编程交流平台!开放,进步,分享!让技术改变生活,让极客改变未来! Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...