# Load Data
path = "X://scJulia/data/10X_PBMC/filtered_feature_bc_matrix/"
TestMat = read_10X("S://PissaroData/10X_PBMC/filtered_feature_bc_matrix/");
B = deepcopy(TestMat["Gene Expression"]);

@time NormalizeData(B.array)

@time pissaro = init_pissaro(B, min_features = 0, min_cells = 5);
@time NormalizeData!(pissaro); #0.17s
@time lul = ScaleData(pissaro, assay = "RNA") #0.6
@time ScaleData!(pissaro, assay = "RNA") #0.6
dimnames(pissaro.assays["RNA"].norm_count)
intersect(names(pissaro.assays["RNA"].norm_count)[1], ["BPREB3", "CD79A", "CD79B", "DNTT"])
peepee = names(lul.assays["RNA"].scaled_count)
pissaro.assays["RNA"].norm_count[["CD79A", "CD79B", "DNTT"],:]
lul.assays["RNA"].scaled_count
pissaro.assays["RNA"].scaled_count
scale_factor = 10000
@time lul, dick, shit = tsvd(pissaro.assays["RNA"].scaled_count, 50)
ScaleData(pissaro.assays["RNA"].norm_count.array)
Base.summarysize(flunch)/(1024^2)
GetNormCount(flunch)
@time scaled_count(flunch)
FindVariableFeatures(flunch)
GetVariableFeatures(flunch)
GetFeatureData(flunch)
PercentageFeatureSet(pissaro)
pissaro.assays["RNA"].scaled_count["CD79A",:]
flunch.assays["RNA"].norm_count.array
A = convert(SparseMatrixCSC,B)
k = 10
m, n = size(A)
F = Arpack.svds(A, nsv = k)[1]
U::Array{Float32, 2}, Z::Array{Float32, 1}, V::Array{Float32, 2} = F.U, F.S, F.V
U, Z, V = convert(Array{Float32, 2}, F.U), convert(Array{Float32, 1}, F.S), convert(Array{Float32, 2}, F.V)
@assert size(U, 2) == length(Z) == size(V, 2) == k
return U, Z, V

B.array

lol = tsvd(B, 6)
lol[1]
using Plots
plot(lol[1][:,1], lol[1][:,2], seriestype = :scatter)

using MultivariateStats
MultivariateStats.ica
ok = fit(ICA, A, 5)

ScaledData = (C .- mean(C,dims=2)) ./ std(C,dims=2)

A1 = [[1]]
A2 = copy(A1)
A1[1] === A2[1]                 # true
A3 = deepcopy(A1)
A1[1] === A3[1]  

C = deepcopy(flunch.assays["RNA"].norm_count.array)
StatsBase.standardize(ZScoreTransform, C.array[1:100, :], dims = 2)

GetCount(flunch, assay = "DICK")
in("DICK",keys(flunch.assays))

function test1(x::Int; y)
    x+y
end

Matrix(scaled_count(flunch))
scale_factor = 10000
GetCount(X, assay = assay)
GetNormCount(X, assay = assay)
NormalizeData(GetCount(flunch, assay = "RNA"), scale_factor = 2000)
X = flunch.assays["RNA"].count
log1p.(scale_factor * (X ./ sum(X, dims = 1)))
GetNormCount(flunch, assay = "RNA", named = false)

@time for i in ProgressBar(1:10000000) #wrap any iterator
    sqrt(i)
 end
@time sqrt.(1:(10^9))

function copy_cols(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{T}, inds, inds)
    for i = inds
        out[:, i] = x
    end
    return out
end

function copy_rows(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{T}, inds, inds)
    for i = inds
        out[i, :] = x
    end
    return out
end

function copy_col_row(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{Float64}, inds, inds)
    for col = inds, row = inds
        out[row, col] = x[row]
    end
    return out
end

function copy_row_col(x::Vector{T}) where T
    inds = axes(x, 1)
    out = similar(Array{T}, inds, inds)
    for row = inds, col = inds
        out[row, col] = x[col]
    end
    return out
end


x = randn(10000);

@time copy_cols(x)';
@time copy_col_row(x);
@time copy_rows(x);

fmt(f) = println(rpad(string(f)*": ", 14, ' '), @elapsed f(x))

map(fmt, [copy_cols, copy_rows, copy_col_row, copy_row_col]);


df = DataFrame(val=rand(8), injury = (repeat(["shoulder", "leg"],4)))


injuries["shoulder"]
A = feat
B = unique(A)
for i in B
    A[A .== i] = A[A .== i] .* append!([""],(string.(range(1,stop =length(A[A .== i])-1, step = 1)).* "."))
end

dict = Dict(i for i in A)
A = ["dog", "dog", "cat", "dog", "cat"]

function make_unique(x::Array{String,1})::Array{String,1}
    dict = Dict{String, Int}()
    for j in unique(x)
        dict[j] = 0
    end    
    y = Array{Int,1}()
    for i in x
        push!(y, dict[i])
        dict[i] = dict[i]+1
    end
    x .*  replace("."  .* string.(y), ".0" => "")
end

mmfile = GZip.open(rmat_path)
@time s = deepcopy(readlines(fh));

open(path * "bite.mtx", "w") do io
    writedlm(io, readlines(fh))
end

open(path * "bite.mtx","w") do io
    write(io,s)
end

writedlm(path * "bite.mtx", readlines(fh))

@time lol = mmread(rmat_path);
@time lol2 = mmread2(rmat_path .* ".gz")
mmfile = GZip.open(rmat_path, "r")


function make_unique(x::Array{String,1}; ram::Bool = true)
    if(ram)
	    r=similar(x)
	    count=1
	    sp=sortperm(x)
	    r[sp[1]]=x[sp[1]]
	    for index in 2:length(x)
	    	if x[sp[index]]==x[sp[index-1]]
	    		r[sp[index]]=x[sp[index]]*"."*string(count)
	    		count+=1
	    	else
	    		r[sp[index]]=x[sp[index]]
	    		count=1
	    	end
	    end
	    return r
    else    
        counts = Dict{String, Int}()
    result = Vector{String}()
    for e in x
        if e in keys(counts)
            push!(result, string(e, ".", counts[e]))
            counts[e] += 1
        else
            push!(result, e)
            counts[e] = 1
        end
    end
    return result
    end
end

pissaro.assays["RNA"].count.array = 
assays["RNA"].scaled_count

span = 0.3

function FindVariableFeatures(mat::AbstractArray; selection_method::String = "vst", loess_span::Float64 = 0.3, 
    clip_max::Union{String, Float64} = "auto", mean_function::Function = FastExpMean, dispersion_function::Function = FastLogVMR, 
    num_bin::Int64 = 20, binning_method::String = "equal_width")::DataFrame
    mat = pissaro.assays["RNA"].count.array
    if clip_max == "auto"
            clip_max = sqrt(mat.n)
    end
    hvf_mean = mean(mat, dims = 2)[:,1]
    hvf_var  = var(mat, dims = 2)[:,1]
    hvf_expected = Vector{Float64}(undef, length(hvf_var))
    hvf_not_const = hvf_var .> 0 
    res = loess(log10.(hvf_mean[hvf_not_const]), log10.(hvf_var[hvf_not_const]), span = 0.3)
    hvf_expected[hvf_not_const] .= 10 .^(Loess.predict(res, log10.(hvf_mean)))
    a, b, c = findnz(mat)
    nonzero_var =  min.((c .- hvf_mean[a]) ./ sqrt.(hvf_expected)[a], clip_max) .^2
    df = DataFrame(A = a, B = b, C = c, D = nonzero_var)
    df = sort(combine(DataFrames.groupby(df, :A), :D => sum), :A)
    hvf_std = (df.D_sum .+ zero_var)/(mat.n - 1)
    DataFrame(vst_mean = hvf_mean, vst_variance = hvf_var, vst_variance_expected = hvf_expected, vst_variance_standardized = hvf_std)
end


@time lul =  FindVariableFeatures(pissaro)
@time FindVariableFeatures!(pissaro)
FindVariableFeatures(B.assays["RNA"].count)
fieldnames(typeof(lul.assays["RNA"]))  
(sort(pissaro.assays["RNA"].feature_data, :vst_variance_standardized, rev = true).features)[1:2000]
B = deepcopy(pissaro)
DataFrames.select(B.assays["RNA"].feature_data, DataFrames.Not(["vst_mean", "vst_variance", "vst_variance_expected", "vst_variance_standardized"]))
var_col = ["vst_mean", "vst_variance", "vst_variance_expected", "vst_variance_standardized"]

function our_type(number)
    add5() = number += 5
    double() = number *= 2
    get_number() = number
    (var) -> (add5;double;number; get_number)
end

bite = our_type(10)
bite.add5()
bite.get_number()