Basic X0 model

documentaiton of this will be added soon

library(GenomeTelescope)

kmer_hist_female <- read.table('data/springtails/WW5-3_k21_truncated.hist', header = T)

female_model <- genomescope(kmer_hist_female, 'data/springtails/WW5-3')

female_coefficients_est <- coef(female_model[[1]])

kmer_hist_male <- read.table('data/springtails/Afus1_k21_truncated.hist', header = T)

cov_range <- 15:150
x <- kmer_hist_male[cov_range, 1]
y <- kmer_hist_male[cov_range, 2]

naive_kmer_est <- x[which.max(y)] / 2
lengthEst <- sum(x * y) / naive_kmer_est
hetEst <- 0.6

two_tissue_model <- nlsLM_2peak_two_tissue_model(x, y, naive_kmer_est, lengthEst, hetEst)

plot_two_tissue_model(two_tissue_model)
two_tissue_est <- coef(two_tissue_model)

sperm_coverage <- coef(two_tissue_model)['kmercov'] - (coef(two_tissue_model)['kmercov2'] - coef(two_tissue_model)['kmercov'])
x <- round(x - sperm_coverage)

estKmercov <- coef(two_tissue_model)['kmercov2'] - coef(two_tissue_model)['kmercov']
estR <- female_coefficients_est['r']
estLength <- female_coefficients_est['length']

male_model <- nls_2peak_male_X0(x, y, k = 21, estKmercov, estLength, estR)