# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .jl
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Julia 1.6.1
#     language: julia
#     name: julia-1.6
# ---

DATE = "2021-06-24"
TASK = "simplified-kmer-graph"
DIR = "$(DATE)-$(TASK)"
DIR = mkpath("$(homedir())/$(DIR)")

# +
pkgs = [
"Graphs",
"MetaGraphs",
"BioSequences",
"uCSV",
"DataFrames",
"FASTX",
"Random",
"ProgressMeter",
"Revise"
]

import Pkg
Pkg.add(pkgs)
for pkg in pkgs
    eval(Meta.parse("import $(basename(pkg))"))
end

import Mycelia
# -

# set a random seed
seed = Random.seed!(0)

# randomly generate a dna sequence of 100bp
genome = BioSequences.randdnaseq(seed, 100)

# define error rate
error_rate = 0.01

# generate 100x coverage fastq file
# put accuracy rate into fastq file
coverage = 10
fastq_file = "$(DIR)/$(DATE)-$(TASK).fastq"
open(fastq_file, "w") do io
    fastq_writer = FASTX.FASTQ.Writer(io)
    for i in 1:coverage
        observed_sequence = Mycelia.observe(genome, error_rate=error_rate)
        q = -10 * log10(error_rate)
        quality_scores = fill(q, length(observed_sequence))
        fastq_record = FASTX.FASTQ.Record("i", observed_sequence, quality_scores)
        write(fastq_writer, fastq_record)
    end
end

k = 11

kmer_type = BioSequences.BigDNAMer{k}

simple_kmer_graph = Mycelia.fastx_to_simple_kmer_graph(kmer_type, fastq_file)

# visualize
gfa_file = fastq_file * ".k-$k.gfa"
Mycelia.graph_to_gfa(simple_kmer_graph, gfa_file)

run(`/Applications/Bandage.app/Contents/MacOS/Bandage image $(gfa_file) $(gfa_file).svg --depwidth 1 --deppower 1`)
# --nodewidth <float> Average node width (0.5 to 1000, default: 5)
# --depwidth <float>  Depth effect on width (0 to 1, default: 0.5)
# --deppower <float>  Power of depth effect on width (0 to 1, default: 0.5)

html_path_to_svg = "./" * repeat("../", length(split(pwd(), '/')) - 3)
html_path_to_svg *= replace("$(gfa_file).svg", "$(homedir())/" => "")

x = display("text/html", "<img src=$(html_path_to_svg)>")

output_fastq_file = Mycelia.polish_fastq(simple_kmer_graph, fastq_file)

simple_kmer_graph = Mycelia.fastx_to_simple_kmer_graph(kmer_type, output_fastq_file)

# visualize
gfa_file = output_fastq_file * ".k-$k.gfa"
Mycelia.graph_to_gfa(simple_kmer_graph, gfa_file)

run(`/Applications/Bandage.app/Contents/MacOS/Bandage image $(gfa_file) $(gfa_file).svg --depwidth 1 --deppower 1`)
# --nodewidth <float> Average node width (0.5 to 1000, default: 5)
# --depwidth <float>  Depth effect on width (0 to 1, default: 0.5)
# --deppower <float>  Power of depth effect on width (0 to 1, default: 0.5)

html_path_to_svg = "./" * repeat("../", length(split(pwd(), '/')) - 3)
html_path_to_svg *= replace("$(gfa_file).svg", "$(homedir())/" => "")

x = display("text/html", "<img src=$(html_path_to_svg)>")

StarcoderdataJupyterScriptsDedupFiltered