Correlation_plot

ggplotRegression <- function (title, fit) {

require(ggplot2)
print(summary(fit))
ggplot(fit$model, aes_string(x = names(fit$model)[2], y = names(fit$model)[1])) + 
  geom_point() +
  stat_smooth(method = "lm", col = "red") +
  labs(title = paste(title, " Adj R2 = ",signif(summary(fit)$adj.r.squared, 5)), x = "rep1 log2 reads", y = "rep2 log2 reads")
}

10k correlation

rep1 <- read.table("10k_leaf_ATAC_rep1_readcount.txt")
colnames(rep1) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")
rep2 <- read.table("10k_leaf_ATAC_rep2_readcount.txt")
colnames(rep2) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")

reps <- data.frame(rep1 = log2(rep1$read_count), rep2 = log2(rep2$read_count))
fit <- lm(rep2 ~ rep1, data = reps)
ggplotRegression("10k", fit)

## Loading required package: ggplot2

## 
## Call:
## lm(formula = rep2 ~ rep1, data = reps)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.77874 -0.18858  0.00936  0.19816  2.33138 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.756083   0.011798   148.8   <2e-16 ***
## rep1        1.011328   0.001951   518.3   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.32 on 14167 degrees of freedom
## Multiple R-squared:  0.9499, Adjusted R-squared:  0.9499 
## F-statistic: 2.686e+05 on 1 and 14167 DF,  p-value: < 2.2e-16

## `geom_smooth()` using formula 'y ~ x'

20k correlation

rep1 <- read.table("20k_leaf_ATAC_rep1_readcount.txt")
colnames(rep1) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")
rep2 <- read.table("20k_leaf_ATAC_rep2_readcount.txt")
colnames(rep2) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")

reps <- data.frame(rep1 = log2(rep1$read_count), rep2 = log2(rep2$read_count))
fit <- lm(rep2 ~ rep1, data = reps)
ggplotRegression("20k", fit)

## 
## Call:
## lm(formula = rep2 ~ rep1, data = reps)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.81868 -0.19948  0.00827  0.20774  1.48079 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.472723   0.019688   -74.8   <2e-16 ***
## rep1         0.887915   0.002453   362.0   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3238 on 9013 degrees of freedom
## Multiple R-squared:  0.9356, Adjusted R-squared:  0.9356 
## F-statistic: 1.31e+05 on 1 and 9013 DF,  p-value: < 2.2e-16

## `geom_smooth()` using formula 'y ~ x'

50k correlation

rep1 <- read.table("50k_leaf_ATAC_rep1_readcount.txt")
colnames(rep1) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")
rep2 <- read.table("50k_leaf_ATAC_rep2_readcount.txt")
colnames(rep2) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")

reps <- data.frame(rep1 = log2(rep1$read_count), rep2 = log2(rep2$read_count))
fit <- lm(rep2 ~ rep1, data = reps)
ggplotRegression("50k", fit)

## 
## Call:
## lm(formula = rep2 ~ rep1, data = reps)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.16374 -0.18653  0.00508  0.18099  1.62245 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.765484   0.015963  -110.6   <2e-16 ***
## rep1         0.953442   0.001991   478.8   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2995 on 12170 degrees of freedom
## Multiple R-squared:  0.9496, Adjusted R-squared:  0.9496 
## F-statistic: 2.293e+05 on 1 and 12170 DF,  p-value: < 2.2e-16

## `geom_smooth()` using formula 'y ~ x'

80k correlation

rep1 <- read.table("80k_leaf_ATAC_rep1_readcount.txt")
colnames(rep1) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")
rep2 <- read.table("80k_leaf_ATAC_rep2_readcount.txt")
colnames(rep2) <- c("chr", "start", "end", "read_count", "coverage", "length", "fraction")

reps <- data.frame(rep1 = log2(rep1$read_count), rep2 = log2(rep2$read_count))
fit <- lm(rep2 ~ rep1, data = reps)
ggplotRegression("80k", fit)

## 
## Call:
## lm(formula = rep2 ~ rep1, data = reps)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.87473 -0.18865  0.00958  0.19216  1.85952 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.462717   0.016507  -88.61   <2e-16 ***
## rep1         0.922775   0.002084  442.86   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3044 on 10939 degrees of freedom
## Multiple R-squared:  0.9472, Adjusted R-squared:  0.9472 
## F-statistic: 1.961e+05 on 1 and 10939 DF,  p-value: < 2.2e-16

## `geom_smooth()` using formula 'y ~ x'