Skip to content

Commit

Permalink
Format number of poly sites
Browse files Browse the repository at this point in the history
  • Loading branch information
SeviJordi committed Oct 10, 2023
1 parent d943931 commit c1f49ae
Show file tree
Hide file tree
Showing 2 changed files with 10 additions and 6 deletions.
2 changes: 1 addition & 1 deletion template.qmd
Original file line number Diff line number Diff line change
Expand Up @@ -143,7 +143,7 @@ The nucleotide diversity value for the studied samples is $`r div_values[["diver

### Number of polymorphic sites

Sites with minor allele frequency > 0.05 are considered as polymorphic sites.
Sites with minor allele frequency > 0.05 are considered as polymorphic sites. Collection date explains the number of polymorphic sites with an $R^2$ of **`r nv.counts[["r2"]]`** with a p-value of $`r nv.counts[["value"]]`$.

![**Figure 6.** Number of polymorphic sites depending on collection date. The red line shows the linear model fit.](`r params$fig_cor_snp`)

Expand Down
14 changes: 9 additions & 5 deletions workflow/scripts/report/NV_description.R
Original file line number Diff line number Diff line change
Expand Up @@ -360,7 +360,7 @@ ggsave(

# Figure for nº of heterozygus sites for each sample
log_info("Plotting nº of heterozygus sites for each sample")
figur_SNP_time <- vcf_snp %>%
figur_SNP_table <- vcf_snp %>%
filter(ALT_FREQ <= 0.95) %>%
select(!GFF_FEATURE) %>%
unique() %>%
Expand All @@ -373,7 +373,9 @@ figur_SNP_time <- vcf_snp %>%
CollectionDate = min(as.Date(CollectionDate)),
n = n_distinct(POS)
) %>%
ungroup() %>%
ungroup()

figur_SNP_time <- figur_SNP_table %>%
ggplot() +
aes(
x = CollectionDate,
Expand All @@ -385,7 +387,6 @@ figur_SNP_time <- vcf_snp %>%
alpha = 0.6
) +
geom_point() +
stat_cor(geom = "label") +
labs(
x = "Date",
y = "Nº of polimorphic sites"
Expand Down Expand Up @@ -467,12 +468,15 @@ n_indels <- vcf %>%
length()

n_snv <- length(unique(vcf$SNP)) - n_indels

model <- lm(n ~ CollectionDate, data = figur_SNP_table)
cortest <- cor.test(figur_SNP_table$n, as.numeric(figur_SNP_table$CollectionDate))
list(
"INDELS" = n_indels,
"SNV" = n_snv,
"window" = snakemake@params[["window"]],
"step" = snakemake@params[["step"]]
"step" = snakemake@params[["step"]],
"r2" = summary(model)$r.squared[[1]],
"value" = ifelse(cortest$p.value < 0.001, "< 0.001", cortest$p.value)
) %>%
toJSON() %>%
write(snakemake@output[["json"]])

0 comments on commit c1f49ae

Please sign in to comment.