text:
# syntax: proto3
using ProtoBuf
import ProtoBuf.meta

mutable struct TensorShapeProto_Dim <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function TensorShapeProto_Dim(; kwargs...)
        obj = new(meta(TensorShapeProto_Dim), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct TensorShapeProto_Dim
const __meta_TensorShapeProto_Dim = Ref{ProtoMeta}()
function meta(::Type{TensorShapeProto_Dim})
    ProtoBuf.metalock() do
        if !isassigned(__meta_TensorShapeProto_Dim)
            __meta_TensorShapeProto_Dim[] = target = ProtoMeta(TensorShapeProto_Dim)
            allflds = Pair{Symbol,Union{Type,String}}[:size => Int64, :name => AbstractString]
            meta(target, TensorShapeProto_Dim, allflds, ProtoBuf.DEF_REQ, ProtoBuf.DEF_FNUM, ProtoBuf.DEF_VAL, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_TensorShapeProto_Dim[]
    end
end
function Base.getproperty(obj::TensorShapeProto_Dim, name::Symbol)
    if name === :size
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :name
        return (obj.__protobuf_jl_internal_values[name])::AbstractString
    else
        getfield(obj, name)
    end
end

mutable struct TensorShapeProto <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function TensorShapeProto(; kwargs...)
        obj = new(meta(TensorShapeProto), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct TensorShapeProto
const __meta_TensorShapeProto = Ref{ProtoMeta}()
function meta(::Type{TensorShapeProto})
    ProtoBuf.metalock() do
        if !isassigned(__meta_TensorShapeProto)
            __meta_TensorShapeProto[] = target = ProtoMeta(TensorShapeProto)
            fnum = Int[2,3]
            allflds = Pair{Symbol,Union{Type,String}}[:dim => Base.Vector{TensorShapeProto_Dim}, :unknown_rank => Bool]
            meta(target, TensorShapeProto, allflds, ProtoBuf.DEF_REQ, fnum, ProtoBuf.DEF_VAL, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_TensorShapeProto[]
    end
end
function Base.getproperty(obj::TensorShapeProto, name::Symbol)
    if name === :dim
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{TensorShapeProto_Dim}
    elseif name === :unknown_rank
        return (obj.__protobuf_jl_internal_values[name])::Bool
    else
        getfield(obj, name)
    end
end

export TensorShapeProto_Dim, TensorShapeProto


text:
# Ranking based Metrics for Recommender Systems

"""
    ranking_stats_k(y_rec, y_rel, k = 10)

Evaluates the relevancy of top k recommendations using `precison@k`, `recall@k` and `f1_score@k`. Returns result as a `Dict`.

Here, `y_rec` are predicted probabilities for recommendation and `y_rel` defines as `1` if particular result is relevant, else `0`.
The shape of `y_rec` and `y_rel` are expected to be (1, N_elements).<br>

* `precison_k` is evaluated as `Recommended_items_that_are_relevant / Total_Recommended_items`.
* `recall_l` is evaluated as `Recommended_items_that_are_relevant / Total_Relevant_items`.
* `f1_k` is evaluated as `2 * Recommended_items_that_are_relevant / (Total_Recommended_items + Total_Relevant_items)`.
"""
function ranking_stats_k(y_rec, y_rel; k = 10)
    @assert size(y_rec) == size(y_rel)
    if size(y_rec,1) !=1
       y_rec = reshape(y_rec, 1, :)
       y_rel = reshape(y_rel, 1, :)
    end
    top_k = partialsortperm(y_rec[1,:], 1:k, lt= >)
    y_rec_k, y_rel_k = y_rec[:, top_k], y_rel[:, top_k]
    tp = sum(y_rel_k .== 1)
    total_positive = sum(y_rel.==1)
    precision_k = k!=0 ? tp / k : 1
    recall_k = total_positive!=0 ? tp / total_positive : 1
    f1_k = (total_positive + k!=0) ? (2 * tp / (total_positive + k)) : 1
    return Dict("precision_k" => precision_k,
                "recall_k" => recall_k,
                "f1_k" => f1_k)
end

"""
    avg_precision(y_rec, y_rel, k = 10)

Evaluates how much of the relevant documents are concentrated in the highest ranked predictions.

Calculated as
    `∑(Recall@i - Recall@i-1)* Precision@i for i = (1, 2, 3....k)`

Here, `y_rec` are predicted probabilities for recommendation and `y_rel` defines as `1` if particular result is relevant, else `0`.
The shape of `y_rec` and `y_rel` are expected to be (1, N_elements)
"""
function avg_precision(y_rec, y_rel; k = 10)
    @assert size(y_rec) == size(y_rel)
    if size(y_rec,1) !=1
       y_rec = reshape(y_rec, 1, :)
       y_rel = reshape(y_rel, 1, :)
    end
    top_k = partialsortperm(y_rec[1,:], 1:k, lt= >)
    y_rec_k, y_rel_k = y_rec[:, top_k], y_rel[:, top_k]
    correct_prediction = 0
    running_sum = 0
    for i in 1:k
        if y_rel_k[i] == 1
            correct_prediction +=1
            running_sum += correct_prediction / i
        end
    end
    return running_sum / sum(y_rel .== 1)
end

# TODO

# DCG
# NDCG

ext: jl


text:
using Test
using MPI

if get(ENV,"JULIA_MPI_TEST_ARRAYTYPE","") == "CuArray"
    import CUDA
    ArrayType = CUDA.CuArray
else
    ArrayType = Array
end

MPI.Init()

comm = MPI.COMM_WORLD
size = MPI.Comm_size(comm)
rank = MPI.Comm_rank(comm)

dst = mod(rank+1, size)
src = mod(rank-1, size)

N = 32

send_mesg = ArrayType{Float64}(undef,N)
recv_mesg = ArrayType{Float64}(undef,N)
recv_mesg_expected = ArrayType{Float64}(undef,N)

fill!(send_mesg, Float64(rank))
fill!(recv_mesg_expected, Float64(src))

rreq = MPI.Irecv!(recv_mesg, comm; source=src,  tag=src+32)
sreq = MPI.Isend(send_mesg, comm; dest=dst, tag=rank+32)

stats = MPI.Waitall([sreq, rreq], MPI.Status)
@test rreq isa MPI.Request
@test sreq isa MPI.Request
@test MPI.Get_source(stats[2]) == src
@test MPI.Get_tag(stats[2]) == src+32
@test recv_mesg == recv_mesg_expected

@test MPI.Testall([sreq, rreq])

if rank == 0
    MPI.send(send_mesg, comm; dest=dst, tag=rank+32)
    recv_mesg = recv_mesg_expected
elseif rank == size-1
    recv_mesg = MPI.recv(comm; source=src, tag=src+32)
else
    recv_mesg = MPI.recv(comm; source=src, tag=src+32)
    MPI.send(send_mesg, comm; dest=dst, tag=rank+32)
end
@test recv_mesg == recv_mesg_expected

if rank == 0
    MPI.Send(Float64(rank), comm; dest=dst, tag=rank+32)
    recv_val = Float64(src)
elseif rank == size-1
    recv_val = MPI.Recv(Float64, comm; source=src, tag=src+32)
else
    recv_val = MPI.Recv(Float64, comm; source=src, tag=src+32)
    MPI.Send(Float64(rank), comm; dest=dst, tag=rank+32)
end
@test recv_val == Float64(src)

rreq = MPI.Irecv!(recv_mesg, comm; source=src,  tag=src+32)
sreq = MPI.Isend(send_mesg, comm; dest=dst, tag=rank+32)

req_arr = [sreq,rreq]
inds = MPI.Waitsome(req_arr)
for i in inds
    @test MPI.Test(req_arr[i])
end

rreq = MPI.Irecv!(recv_mesg, comm; source=src,  tag=src+32)
MPI.Cancel!(rreq)
MPI.Wait(rreq)
@test rreq.buffer == nothing

GC.gc()

# ---------------------
# MPI_Sendrecv function
# ---------------------
#
# send datatype
# ---------------------
# We test this function by executing a left shift of the leftmost element in a 1D
# cartesian topology with periodic boundary conditions.
#
# Assume we have two processors, the data will look like this
# proc 0 | proc 1
#  0 0 0 |  1 1 1
#
# After the shift the data will contain
# proc 0 | proc 1
#  0 0 1 |  1 1 0
#
# init data
comm_rank = MPI.Comm_rank(comm)
comm_size = MPI.Comm_size(comm)
a = Float64[comm_rank, comm_rank, comm_rank]

# construct cartesian communicator with 1D topology
comm_cart = MPI.Cart_create(comm, 1, Cint[comm_size], Cint[1], false)

# get source and dest ranks using Cart_shift
src_rank, dest_rank = MPI.Cart_shift(comm_cart, 0, -1)

# execute left shift using subarrays
MPI.Sendrecv!(@view(a[1]), @view(a[3]), comm_cart;
        dest=dest_rank, sendtag=0, source=src_rank, recvtag=0)

@test a == [comm_rank, comm_rank, (comm_rank+1) % comm_size]

# send elements from a buffer
# ---------------------------
a = Float64[comm_rank, comm_rank, comm_rank]
b = Float64[       -1,        -1,        -1]
MPI.Sendrecv!(@view(a[1:2]), @view(b[1:2]), comm_cart;
        dest=dest_rank, sendtag=1, source=src_rank, recvtag=1)
@test b == [(comm_rank+1) % comm_size, (comm_rank+1) % comm_size, -1]

# send entire buffer
# ---------------------------
a = Float64[comm_rank, comm_rank, comm_rank]
b = Float64[       -1,        -1,        -1]
MPI.Sendrecv!(a, b,  comm_cart;
        dest=dest_rank, sendtag=2, source=src_rank, recvtag=2)

@test b == [(comm_rank+1) % comm_size, (comm_rank+1) % comm_size, (comm_rank+1) % comm_size]

@test MPI.Waitall(MPI.Request[], MPI.Status) == MPI.Status[]

MPI.Finalize()
# @test MPI.Finalized()


text:
using Random
export RandomShuffle

"""
    RandomShuffle() <: Surrogate

A random constrained surrogate, generated by shifting values around.

Random shuffle surrogates preserve the mean, variance and amplitude
distribution of the original signal. Properties not preserved are *any
temporal information*, such as the power spectrum and hence linear
correlations.

The null hypothesis this method can test for is whether the data
are uncorrelated noise, possibly measured via a nonlinear function.
Specifically, random shuffle surrogate can test
the null hypothesis that the original signal is produced by independent and
identically distributed random variables[^Theiler1991, ^Lancaster2018].

*Beware: random shuffle surrogates do not cover the case of correlated noise*[^Lancaster2018].

[^Theiler1991]: J. Theiler, S. Eubank, A. Longtin, B. Galdrikian, J. Farmer, Testing for nonlinearity in time series: The method of surrogate data, Physica D 58 (1–4) (1992) 77–94.
"""
struct RandomShuffle <: Surrogate end

function surrogenerator(x::AbstractVector, rf::RandomShuffle, rng = Random.default_rng())
    n = length(x)
    idxs = collect(1:n)

    init = (
        permutation = zeros(Int, n),
        idxs = idxs,
    )

    return SurrogateGenerator(rf, x, similar(x), init, rng)
end

function (sg::SurrogateGenerator{<:RandomShuffle})()
    x, s, rng = sg.x, sg.s, sg.rng
    permutation, idxs  = getfield.(Ref(sg.init), (:permutation, :idxs))
    sample!(rng, idxs, permutation; replace = false)
    s .= x[permutation]
    return s
end


text:
# Analyse results

export
    writepdb

# Extension of method from BioStructures
function BioStructures.writepdb(filepath::AbstractString, universe::Universe)
    open(filepath, "w") do output
        for (i, c) in enumerate(universe.coords)
            at = universe.molecule.atoms[i]
            at_rec = AtomRecord(
                false,
                i,
                at.name,
                ' ',
                at.resname,
                "A",
                at.resnum,
                ' ',
                [10*c.x, 10*c.y, 10*c.z],
                1.0,
                0.0,
                "  ",
                "  "
            )
            println(output, pdbline(at_rec))
        end
    end
end


text:
module TestLookup
using CompressHashDisplace
using Test

@testset "FrozenDict Lookup" begin
    dict = Dict{String, Int}("hello" => 1, "world" => 2)
    fd = FrozenDict(dict)
    @test fd["hello"] == 1
    @test fd["world"] == 2
    @test_throws KeyError fd["goodbye"]

    fud = FrozenUnsafeDict(dict)
    @test fd["hello"] == 1
    @test fd["world"] == 2
end

@testset "Nonstring Lookup" begin
    dict = Dict{Int, Int}(5 => 1, 6 => 2)
    fd = FrozenDict(dict)
    @test fd[5] == 1
    @test fd[6] == 2
    @test_throws KeyError fd["goodbye"]

    fud = FrozenUnsafeDict(dict)
    @test fd[5] == 1
    @test fd[6] == 2
end

end # module


text:
using ..DGmethods.NumericalFluxes:
    NumericalFluxNonDiffusive, NumericalFluxGradient, NumericalFluxDiffusive

import ..DGmethods: boundary_state!

"""
    boundary_state!(nf, ::HBModel, args...)

applies boundary conditions for the hyperbolic fluxes
dispatches to a function in OceanBoundaryConditions.jl based on bytype defined by a problem such as SimpleBoxProblem.jl
"""
@inline function boundary_state!(nf, ocean::HBModel, args...)
    boundary_condition = ocean.problem.boundary_condition
    return ocean_boundary_state!(nf, boundary_condition, ocean, args...)
end

"""
    boundary_state!(nf, ::LinearHBModel, args...)

applies boundary conditions for the hyperbolic fluxes
dispatches to a function in OceanBoundaryConditions.jl based on bytype defined by a problem such as SimpleBoxProblem.jl
"""
@inline function boundary_state!(nf, linear::LinearHBModel, args...)
    ocean = linear.ocean
    boundary_condition = ocean.problem.boundary_condition

    return ocean_boundary_state!(nf, boundary_condition, ocean, args...)
end

"""
    ocean_boundary_state!(nf, boundaries::Tuple, ::HBModel,
                          Q⁺, A⁺, n, Q⁻, A⁻, bctype)
applies boundary conditions for the first-order and gradient fluxes
dispatches to a function in OceanBoundaryConditions.jl based on bytype defined by a problem such as SimpleBoxProblem.jl
"""
@generated function ocean_boundary_state!(
    nf::Union{NumericalFluxNonDiffusive, NumericalFluxGradient},
    boundaries::Tuple,
    ocean,
    Q⁺,
    A⁺,
    n,
    Q⁻,
    A⁻,
    bctype,
    t,
    args...,
)
    N = fieldcount(boundaries)
    return quote
        Base.Cartesian.@nif(
            $(N + 1),
            i -> bctype == i, # conditionexpr
            i -> ocean_boundary_state!(
                nf,
                boundaries[i],
                ocean,
                Q⁺,
                A⁺,
                n,
                Q⁻,
                A⁻,
                t,
            ), # expr
            i -> error("Invalid boundary tag")
        ) # elseexpr
        return nothing
    end
end

"""
    ocean_boundary_state!(nf, boundaries::Tuple, ::HBModel,
                          Q⁺, A⁺, D⁺, n, Q⁻, A⁻, D⁻, bctype)
applies boundary conditions for the second-order fluxes
dispatches to a function in OceanBoundaryConditions.jl based on bytype defined by a problem such as SimpleBoxProblem.jl
"""
@generated function ocean_boundary_state!(
    nf::NumericalFluxDiffusive,
    boundaries::Tuple,
    ocean,
    Q⁺,
    D⁺,
    A⁺,
    n,
    Q⁻,
    D⁻,
    A⁻,
    bctype,
    t,
    args...,
)
    N = fieldcount(boundaries)
    return quote
        Base.Cartesian.@nif(
            $(N + 1),
            i -> bctype == i, # conditionexpr
            i -> ocean_boundary_state!(
                nf,
                boundaries[i],
                ocean,
                Q⁺,
                D⁺,
                A⁺,
                n,
                Q⁻,
                D⁻,
                A⁻,
                t,
            ), # expr
            i -> error("Invalid boundary tag")
        ) # elseexpr
        return nothing
    end
end

include("OceanBoundaryConditions.jl")

max_stars_repo_name: yairchn/CLIMA

max_issues_repo_name: yairchn/CLIMA

max_forks_repo_name: yairchn/CLIMA


text:
using OffsetArrays

"""
Legendre polynomial on the interval [-1..+1].
"""
@inline function _legendre(j::Int, x::T) where {T}
    if j <= 0
        return one(T)
    elseif j == 1
        return x
    else
        return ( (2j-1) * _legendre(j-1, x) * x - (j-1) * _legendre(j-2, x) ) / j
    end
end

"""
Legendre basis on the interval [0..+1].
"""
struct Legendre{T, LT} <: Basis{T}
    b::LT
    n::Int

    function Legendre{T}(n::Integer) where {T}
        p = n-1
        b = OffsetArray([y -> _legendre(i, 2y-1) * sqrt(T(2i+1)) for i in 0:p], 0:p)
        new{T, typeof(b)}(b, n)
    end

end

Legendre(::Type{T}, n::Integer) where {T} = Legendre{T}(n)
Legendre(n::Integer) = Legendre(Float64, n)

function _eval(L::Legendre{LT}, x::DT, j::Int) where {LT,DT}
    local T = promote_type(LT, DT)
    @assert j ≥ 0 && j < L.n
    _legendre(j, 2x-1) * sqrt(T(2j+1))
end

(L::Legendre)(x::Number, j::Integer) = L.b[j](x)

basis(L::Legendre) = L.b
nbasis(L::Legendre) = L.n
order(L::Legendre) = nbasis(L)
degree(L::Legendre) = nbasis(L) - 1

Base.eltype(::Legendre{T}) where {T} = T
Base.eachindex(L::Legendre) = eachindex(L.b)
Base.axes(L::Legendre) = (Inclusion(0..1), eachindex(L))

Base.hash(L::Legendre, h::UInt) = hash(L.n, h)
Base.:(==)(L1::Legendre, L2::Legendre) = (L1.n == L2.n)
Base.isequal(L1::Legendre{T1}, L2::Legendre{T2}) where {T1,T2} = (T1 == T2 && L1 == L2)

Base.getindex(L::Legendre, x::Number, j::Integer) = L(x,j)
Base.getindex(L::Legendre, x::Number,  ::Colon) = [b(x) for b in L.b]
Base.getindex(L::Legendre, X::AbstractVector, j::Integer) = L.(X,j)
Base.getindex(L::Legendre, X::AbstractVector,  ::Colon) = [b(x) for x in X, b in L.b]

## Derivative

"""
Derivative of Legendre polynomial on the interval [-1..+1].
"""
@inline function _legendre_derivative(j::Int, x::T) where {T}
    if j <= 0
        return zero(T)
    elseif j == 1
        return one(T)
    else
        return ( (2j-1) * _legendre(j-1, x) + (2j-1) * _legendre_derivative(j-1, x) * x - (j-1) * _legendre_derivative(j-2, x) ) / j
    end
end

@simplify *(D::Derivative, L::Legendre) = Mul(D,L)

const LegendreDerivative = QMul2{<:Derivative,<:Legendre}

"""
Evaluate derivative of Legendre polynomial on the interval [0..+1].
"""
function _eval(D::LegendreDerivative, x::DT, j::Int) where {DT}
    @assert j ≥ 0 && j < nbasis(D.B)
    _legendre_derivative(j, promote_type(eltype(D.B), DT)(2x-1)) * 2
end

Base.getindex(D::LegendreDerivative, x::Number, j::Integer) = _eval(D, x, j)
Base.getindex(D::LegendreDerivative, x::Number,  ::Colon) = [_eval(D, x, j) for j in eachindex(D.B)]
Base.getindex(D::LegendreDerivative, X::AbstractVector, j::Integer) = [_eval(D, x, j) for x in X]
Base.getindex(D::LegendreDerivative, X::AbstractVector,  ::Colon) = [_eval(D, x, j) for x in X, j in eachindex(D.B)]

Base.adjoint(L::Legendre) = Derivative(axes(L,1)) * L


text:
@testset "Bathymetry models" begin
    @test true
end


text:
module GarishPrint

export pprint, pprint_struct
using Configurations

# 1.0 Compat
@static if !@isdefined(isnothing)
    isnothing(x) = x === nothing
end

@static if VERSION ≥ v"1.6"
    using Preferences
end

@static if VERSION ≥ v"1.6"
    const color_prefs_toml = @load_preference("color")
else
    const color_prefs_toml = nothing
end

include("color.jl")
include("io.jl")
include("pprint.jl")
include("struct.jl")
include("misc.jl")
include("numbers.jl")
include("arrays.jl")
include("dict.jl")
include("set.jl")

end


text:
_testdata() = artifact"TestData"

cacheddatasets = Dict{String,Arrow.Table}()
"""
    dataset(nm)

Return, as an `Arrow.Table`, the test data set named `nm`, which can be a `String` or `Symbol`
"""
function dataset(nm::AbstractString)
    get!(cacheddatasets, nm) do
        path = joinpath(_testdata(), nm * ".arrow")
        if !isfile(path)
            throw(
                ArgumentError(
                    "Dataset \"$nm\" is not available.\nUse MixedModels.datasets() for available names.",
                ),
            )
        end
        Arrow.Table(path)
    end
end
dataset(nm::Symbol) = dataset(string(nm))

"""
    datasets()

Return a vector of names of the available test data sets
"""
function datasets()
    return first.(
        Base.Filesystem.splitext.(filter(endswith(".arrow"), readdir(_testdata())))
    )
end


text:
struct ImageCreateInfo
    mHandleRef::Ref{vk.VkImageCreateInfo}
    mViewType::vk.VkImageViewType

    function ImageCreateInfo(;
        pNext = C_NULL, #::Ptr{Cvoid}
        flags = 0, #::VkImageCreateFlags
        imageType = vk.VK_IMAGE_TYPE_1D, #::VkImageType
        format = vk.VK_FORMAT_UNDEFINED, #::VkFormat
        extent = vk.VkExtent3D(0, 0, 0), #::VkExtent3D
        mipLevels = 0, #::UInt32
        arrayLayers = 0, #::UInt32
        samples = vk.VK_SAMPLE_COUNT_1_BIT, #::VkSampleCountFlagBits
        tiling = vk.VK_IMAGE_TILING_OPTIMAL, #::VkImageTiling
        usage = vk.VkImageUsageFlags(0), #::VkImageUsageFlags
        sharingMode = vk.VK_SHARING_MODE_EXCLUSIVE, #::VkSharingMode
        queueFamilyIndexCount = 0, #::UInt32
        pQueueFamilyIndices = C_NULL, #::Ptr{UInt32}
        initialLayout = vk.VK_IMAGE_LAYOUT_UNDEFINED, #::VkImageLayout
        imageViewType = vk.VK_IMAGE_VIEW_TYPE_1D #vk.VkImageViewType
    )
        this = new(Ref(vk.VkImageCreateInfo(
            vk.VK_STRUCTURE_TYPE_IMAGE_CREATE_INFO,
            pNext, #Ptr{Cvoid}
            flags, #VkImageCreateFlags
            imageType, #VkImageType
            format, #VkFormat
            extent, #VkExtent3D
            mipLevels, #UInt32
            arrayLayers, #UInt32
            samples, #VkSampleCountFlagBits
            tiling, #VkImageTiling
            usage, #VkImageUsageFlags
            sharingMode, #VkSharingMode
            queueFamilyIndexCount, #UInt32
            pQueueFamilyIndices, #Ptr{UInt32}
            initialLayout #VkImageLayout
        )), imageViewType)
    end
end
@class ImageCreateInfo

function handleRef(this::ImageCreateInfo)::Ref{vk.VkImageCreateInfo}
    return this.mHandleRef
end

function createImageCreateInfo(vkPhyDevice::vk.VkPhysicalDevice;
                               usage = vk.VkImageUsageFlags(0),
                               others...)
    usage |= vk.VK_IMAGE_USAGE_TRANSFER_DST_BIT | vk.VK_IMAGE_USAGE_TRANSFER_SRC_BIT
    tilingRef = Ref{vk.VkImageTiling}(vk.VK_IMAGE_TILING_OPTIMAL)
    if (vkPhyDevice != C_NULL)
        adaptTiling(vkPhyDevice, others[:format], usage, tilingRef)
    end

    return ImageCreateInfo(;initialLayout = vk.VK_IMAGE_LAYOUT_UNDEFINED,
                            mipLevels = 1, arrayLayers = 1,
                            samples = vk.VK_SAMPLE_COUNT_1_BIT,
                            sharingMode = vk.VK_SHARING_MODE_EXCLUSIVE,
                            tiling = tilingRef[],
                            usage = usage,
                            others...)
end

function getCombinedType(viewType::vk.VkImageViewType)
    imgType = 0
    flags = 0
    if viewType == vk.VK_IMAGE_VIEW_TYPE_1D || viewType == vk.VK_IMAGE_VIEW_TYPE_1D_ARRAY
        imgType = vk.VK_IMAGE_TYPE_1D
    elseif viewType == vk.VK_IMAGE_VIEW_TYPE_CUBE || viewType == vk.VK_IMAGE_VIEW_TYPE_CUBE_ARRAY
        imgType = vk.VK_IMAGE_TYPE_2D
        flags |= vk.VK_IMAGE_CREATE_CUBE_COMPATIBLE_BIT
    elseif viewType == vk.VK_IMAGE_VIEW_TYPE_2D || viewType == vk.VK_IMAGE_VIEW_TYPE_2D_ARRAY
        imgType = vk.VK_IMAGE_TYPE_2D
    elseif viewType == vk.VK_IMAGE_VIEW_TYPE_3D
        imgType = vk.VK_IMAGE_TYPE_3D
    end

    return viewType, imgType, flags
end

function adaptTiling(vkPhyDevice::vk.VkPhysicalDevice, format::vk.VkFormat, usage::vk.VkImageUsageFlags, tilingRef::Ref{vk.VkImageTiling})
    features = usageToFeatures(usage)
    props = Ref{vk.VkFormatProperties}()
    vk.vkGetPhysicalDeviceFormatProperties(vkPhyDevice, format, props)
    canUseOptimal = ((props[].optimalTilingFeatures & features) == features)
    canUseLinear = ((props[].linearTilingFeatures & features) == features)

    if (!canUseOptimal && tilingRef[] == vk.VK_IMAGE_TILING_OPTIMAL)
        if (canUseLinear)
            println("The format selected for this image does not support optimal ",
                    "tiling for the required features, falling back to linear ",
                    "tiling.\n",
                    "In order to optimize performance, use a format that supports ",
                    "optimal tiling (e.g. use RGBA instead of RGB textures for ",
                    "attachments)\n",
                    "You can supress this warning by setting the tiling to linear ",
                    "yourself.")
            tilingRef[] = vk.VK_IMAGE_TILING_LINEAR
        else
            error("The format selected can't provide the requested ",
                    "usages. Try another format (e.g. use RGBA ",
                    "instead of RGB for attachments)")
        end
    elseif (!canUseLinear && tilingRef[] == vk.VK_IMAGE_TILING_LINEAR)
        if (canUseOptimal)
            println("The format selected for this image does not support linear ",
                    "tiling for the required features, falling back to optimal ",
                    "tiling.\n",
                    "You can supress this warning by setting the tiling to optimal ",
                    "yourself.")
            tilingRef[] = vk.VK_IMAGE_TILING_OPTIMAL
        else
            error("The format selected can't provide the requested ",
                    "usages. Try another format (e.g. use RGBA ",
                    "instead of RGB for attachments)")
        end
    end
end

function mipLevels(dim::UInt32)::UInt32
    count::UInt32 = UInt32(1)
    while (dim >>= 1) > 0
        count += 1
    end
    return count
end

function texture2D(vkPhyDevice::vk.VkPhysicalDevice,
                         width::UInt32, height::UInt32, format::vk.VkFormat)::ImageCreateInfo
    viewType, imgType, flags = getCombinedType(vk.VK_IMAGE_VIEW_TYPE_2D)
    usage = vk.VK_IMAGE_USAGE_SAMPLED_BIT

    info = createImageCreateInfo(vkPhyDevice, imageViewType = viewType, imageType = imgType, flags = flags,
                                format = format,
                                extent = vk.VkExtent3D(width, height, 1),
                                mipLevels = mipLevels(min(width, height)),
                                usage = usage)
    return info
end

function texture2DArray(vkPhyDevice::vk.VkPhysicalDevice,
                              width::UInt32, height::UInt32,
                             layers::UInt32, format::vk.VkFormat)::ImageCreateInfo
    viewType, imgType, flags = getCombinedType(vk.VK_IMAGE_VIEW_TYPE_2D_ARRAY)
    usage = vk.VK_IMAGE_USAGE_SAMPLED_BIT

    info = createImageCreateInfo(vkPhyDevice, imageViewType = viewType, imageType = imgType, flags = flags,
                             format = format,
                             extent = vk.VkExtent3D(width, height, 1),
                             mipLevels = mipLevels(min(width, height)),
                             arrayLayers = layers,
                             usage = usage)
    return info
end

function storageImage2D(vkPhyDevice::vk.VkPhysicalDevice,
                              width::UInt32, height::UInt32, format::vk.VkFormat)::ImageCreateInfo
    viewType, imgType, flags = getCombinedType(vk.VK_IMAGE_VIEW_TYPE_2D)
    usage = vk.VK_IMAGE_USAGE_STORAGE_BIT

    info = createImageCreateInfo(vkPhyDevice, imageViewType = viewType, imageType = imgType, flags = flags,
                               format = format,
                               extent = vk.VkExtent3D(width, height, 1),
                               usage = usage)
    return info
end

function storageImage2DArray(vkPhyDevice::vk.VkPhysicalDevice,
                                   width::UInt32, height::UInt32,
                                  layers::UInt32, format::vk.VkFormat)::ImageCreateInfo
    viewType, imgType, flags = getCombinedType(vk.VK_IMAGE_VIEW_TYPE_2D_ARRAY)
    usage = vk.VK_IMAGE_USAGE_STORAGE_BIT

    info = createImageCreateInfo(vkPhyDevice, imageViewType = viewType, imageType = imgType, flags = flags,
                               format = format,
                               extent = vk.VkExtent3D(width, height, 1),
                               arrayLayers = layers,
                               usage = usage)
    return info
end

function attachment2D(vkPhyDevice::vk.VkPhysicalDevice,
                    width::Integer, height::Integer, format::vk.VkFormat)::ImageCreateInfo
    viewType, imgType, flags = getCombinedType(vk.VK_IMAGE_VIEW_TYPE_2D)
    usage = vk.VK_IMAGE_USAGE_SAMPLED_BIT
    aspects = aspectsOf(format)
    if (aspects & vk.VK_IMAGE_ASPECT_DEPTH_BIT != 0 ||
        aspects & vk.VK_IMAGE_ASPECT_STENCIL_BIT != 0)
        usage |= vk.VK_IMAGE_USAGE_DEPTH_STENCIL_ATTACHMENT_BIT
    end
    if (aspects & vk.VK_IMAGE_ASPECT_COLOR_BIT != 0)
        usage |= vk.VK_IMAGE_USAGE_COLOR_ATTACHMENT_BIT
    end

    info = createImageCreateInfo(vkPhyDevice, imageViewType = viewType, imageType = imgType, flags = flags,
                                format = format,
                                extent = vk.VkExtent3D(width, height, 1),
                                usage = usage)
    return info
end

function attachment2DArray(vkPhyDevice::vk.VkPhysicalDevice,
                                 width::UInt32, height::UInt32,
                                layers::UInt32, format::vk.VkFormat)::ImageCreateInfo
    viewType, imgType, flags = getCombinedType(vk.VK_IMAGE_VIEW_TYPE_2D_ARRAY)
    usage = vk.VK_IMAGE_USAGE_SAMPLED_BIT
    aspects = aspectsOf(format)
    if (aspects & vk.VK_IMAGE_ASPECT_DEPTH_BIT != 0 ||
        aspects & vk.VK_IMAGE_ASPECT_STENCIL_BIT != 0)
        usage |= vk.VK_IMAGE_USAGE_DEPTH_STENCIL_ATTACHMENT_BIT
    end
    if (aspects & vk.VK_IMAGE_ASPECT_COLOR_BIT != 0)
        usage |= vk.VK_IMAGE_USAGE_COLOR_ATTACHMENT_BIT
    end

    info = createImageCreateInfo(vkPhyDevice, imageViewType = viewType, imageType = imgType, flags = flags,
                                format = format,
                                extent = vk.VkExtent3D(width, height, 1),
                                arrayLayers = layers,
                                usage = usage)
    return info
end


text:
struct RotationMatrix{T}
    m::Matrix{T}
end
function RotationMatrix(phi::T) where {T}
    c = cos(phi)
    s = sin(phi)

    m = Matrix{Float64}(undef, 2, 2)
    m[1,1] = c
    m[2,2] = c
    m[1,2] = -s
    m[2,1] = s

    return RotationMatrix(m)
end

function rotation_matrix(phi)
    # R * v -> v_rotated_by_phi
    c = cos(phi)
    s = sin(phi)

    R = [c -s; s c]
    return R
end

function Base.getindex(rot_mat::RotationMatrix, index)
    return rot_mat.m[index]
end

function Base.getindex(rot_mat::RotationMatrix, index...)
    return rot_mat.m[index...]
end

function rotate(rot_mat, p::T) where {T <: AbstractPoint}
    x = rot_mat[1,1] * p.x + rot_mat[1,2] * p.y
    y = rot_mat[2,1] * p.x + rot_mat[2,2] * p.y
    return T(x, y)
end

function rotate(rot_mat, p::T) where {T <: AbstractArray}
    x = rot_mat[1,1] * p[1] + rot_mat[1,2] * p[1]
    y = rot_mat[2,1] * p[1] + rot_mat[2,2] * p[1]
    return [x, y]
end

max_stars_repo_name: HomoModelicus/julia

max_issues_repo_name: HomoModelicus/julia

max_forks_repo_name: HomoModelicus/julia


text:
module BenchSumTranspose

import Random
using BenchmarkTools
using Transducers

const SUITE = BenchmarkGroup()

function iter_sum(xs, acc = false)
    for x in xs
        acc += x
    end
    return acc
end

function man_sum(xs::AbstractMatrix, acc = false)
    for i in axes(xs, 2), j in axes(xs, 1)
        acc += @inbounds xs[j, i]
    end
    return acc
end

Random.seed!(12345)
# for n in [30, 100]
let n = 30
    s1 = SUITE[string(n)] = BenchmarkGroup()

    xs = randn(n, n)'

    s2 = s1["noinit"] = BenchmarkGroup()
    s2["iter"] = @benchmarkable iter_sum($xs)
    s2["man"] = @benchmarkable man_sum($xs)
    s2["xf"] = @benchmarkable foldxl(+, $xs)

    s2 = s1["withinit"] = BenchmarkGroup()
    s2["iter"] = @benchmarkable iter_sum($xs, 0.0)
    s2["man"] = @benchmarkable man_sum($xs, 0.0)
    s2["xf"] = @benchmarkable foldxl(+, $xs; init = 0.0)
end

end  # module
BenchSumTranspose.SUITE


text:
# Use baremodule to shave off a few KB from the serialized `.ji` file
baremodule cuQuantum_jll
using Base
using Base: UUID
import JLLWrappers

JLLWrappers.@generate_main_file_header("cuQuantum")
JLLWrappers.@generate_main_file("cuQuantum", UUID("b75408ef-6fdf-5d74-b65a-7df000ad18e6"))
end  # module cuQuantum_jll


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

mutable struct ApplicationGatewayFrontendIPConfigurationPropertiesFormat <: SwaggerModel
    privateIPAddress::Any # spec type: Union{ Nothing, String } # spec name: privateIPAddress
    privateIPAllocationMethod::Any # spec type: Union{ Nothing, IPAllocationMethod } # spec name: privateIPAllocationMethod
    subnet::Any # spec type: Union{ Nothing, SubResource } # spec name: subnet
    publicIPAddress::Any # spec type: Union{ Nothing, SubResource } # spec name: publicIPAddress
    privateLinkConfiguration::Any # spec type: Union{ Nothing, SubResource } # spec name: privateLinkConfiguration
    provisioningState::Any # spec type: Union{ Nothing, ProvisioningState } # spec name: provisioningState

    function ApplicationGatewayFrontendIPConfigurationPropertiesFormat(;privateIPAddress=nothing, privateIPAllocationMethod=nothing, subnet=nothing, publicIPAddress=nothing, privateLinkConfiguration=nothing, provisioningState=nothing)
        o = new()
        validate_property(ApplicationGatewayFrontendIPConfigurationPropertiesFormat, Symbol("privateIPAddress"), privateIPAddress)
        setfield!(o, Symbol("privateIPAddress"), privateIPAddress)
        validate_property(ApplicationGatewayFrontendIPConfigurationPropertiesFormat, Symbol("privateIPAllocationMethod"), privateIPAllocationMethod)
        setfield!(o, Symbol("privateIPAllocationMethod"), privateIPAllocationMethod)
        validate_property(ApplicationGatewayFrontendIPConfigurationPropertiesFormat, Symbol("subnet"), subnet)
        setfield!(o, Symbol("subnet"), subnet)
        validate_property(ApplicationGatewayFrontendIPConfigurationPropertiesFormat, Symbol("publicIPAddress"), publicIPAddress)
        setfield!(o, Symbol("publicIPAddress"), publicIPAddress)
        validate_property(ApplicationGatewayFrontendIPConfigurationPropertiesFormat, Symbol("privateLinkConfiguration"), privateLinkConfiguration)
        setfield!(o, Symbol("privateLinkConfiguration"), privateLinkConfiguration)
        validate_property(ApplicationGatewayFrontendIPConfigurationPropertiesFormat, Symbol("provisioningState"), provisioningState)
        setfield!(o, Symbol("provisioningState"), provisioningState)
        o
    end
end # type ApplicationGatewayFrontendIPConfigurationPropertiesFormat

const _property_map_ApplicationGatewayFrontendIPConfigurationPropertiesFormat = Dict{Symbol,Symbol}(Symbol("privateIPAddress")=>Symbol("privateIPAddress"), Symbol("privateIPAllocationMethod")=>Symbol("privateIPAllocationMethod"), Symbol("subnet")=>Symbol("subnet"), Symbol("publicIPAddress")=>Symbol("publicIPAddress"), Symbol("privateLinkConfiguration")=>Symbol("privateLinkConfiguration"), Symbol("provisioningState")=>Symbol("provisioningState"))
const _property_types_ApplicationGatewayFrontendIPConfigurationPropertiesFormat = Dict{Symbol,String}(Symbol("privateIPAddress")=>"String", Symbol("privateIPAllocationMethod")=>"IPAllocationMethod", Symbol("subnet")=>"SubResource", Symbol("publicIPAddress")=>"SubResource", Symbol("privateLinkConfiguration")=>"SubResource", Symbol("provisioningState")=>"ProvisioningState")
Base.propertynames(::Type{ ApplicationGatewayFrontendIPConfigurationPropertiesFormat }) = collect(keys(_property_map_ApplicationGatewayFrontendIPConfigurationPropertiesFormat))
Swagger.property_type(::Type{ ApplicationGatewayFrontendIPConfigurationPropertiesFormat }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_ApplicationGatewayFrontendIPConfigurationPropertiesFormat[name]))}
Swagger.field_name(::Type{ ApplicationGatewayFrontendIPConfigurationPropertiesFormat }, property_name::Symbol) =  _property_map_ApplicationGatewayFrontendIPConfigurationPropertiesFormat[property_name]

function check_required(o::ApplicationGatewayFrontendIPConfigurationPropertiesFormat)
    true
end

function validate_property(::Type{ ApplicationGatewayFrontendIPConfigurationPropertiesFormat }, name::Symbol, val)
end


text:
# This file is a part of GERDADeepLearning.jl, licensed under the MIT License (MIT).

__precompile__(false)

module GERDADeepLearning

include("gerda_files.jl")

include("events.jl")
include("hdf5io.jl")
include("rootio.jl")
include("dl_env.jl")

include("gdl_util.jl")
include("extensions.jl")
include("signal_processing.jl")

include("mx_train_util.jl")
include("conv_autoencoder.jl")
include("dnn_classifier.jl")

include("efficiencies.jl")
include("plotting.jl")

end # module


text:
if !isdefined(Main, :CustomStencil)
    include("../module/CustomStencil.jl")
    using Main.CustomStencil
end
include("misc.jl")
include("cpu_stencils.jl")

## Stencil Definition
coefs = Float64.([0.25; 0.125; 0.06])
dense_stencil = @def_stencil_expression(
            @sum(i,-2,2,
                @sum(j,-2,2,
                    @sum(k,-2,2, c[max(abs(i),abs(j),abs(k))]*D[x+i,y+j,z+k]))))
st_def = CreateStencilDefinition(dense_stencil, coefs)
st_inst = NewStencilInstance(st_def)

## Input Data size Definition
radius = 2
nx = 5
ny = 5
nz = 5
dx = 1<<(nx)
dy = 1<<(ny)
dz = 1<<(nz)

t_steps = 3
## GPU Implementation
original_data = CreateData(dx,dy,dz)
data = PadData(radius, original_data)
dout = ApplyStencil(st_inst, data, t_steps)
gpu_out = Array(dout)

## CPU Implementation
cpu_out =  cpu_dense_stencil(original_data, radius, coefs, t_steps)

## Check Result
mean_er = sum(abs.(cpu_out .- gpu_out))./(reduce(*, size(gpu_out)))
@assert(sum(abs.(cpu_out .- gpu_out))./(reduce(*, size(gpu_out))) < 0.001)
## Visual Check
# using Plots
# shift = -10
# Int(round(dz/2))+shift
# heatmap(cout[:,:,Int(round(dz/2))+shift])
# heatmap(o[:,:,Int(round(dz/2))+shift])


text:
using Random, SnpArrays, DataFrames, GLM
using LinearAlgebra, Test, TraitSimulation
using BenchmarkTools, Statistics
Random.seed!(1234)

function generateSPDmatrix(n)
	A = rand(n)
	m = 0.5 * (A * A')
	PDmat = m + (n * Diagonal(ones(n)))
end

function generateRandomVCM(n::Int64, p::Int64, d::Int64, m::Int64)
	# n-by-p design matrix
	X = randn(n, p)

	# p-by-d mean component regression coefficient for each trait
	B = hcat(ones(p, 1), rand(p))

	V = ntuple(x -> zeros(n, n), m)
	for i = 1:m-1
	  copy!(V[i], generateSPDmatrix(n))
	end
	copy!(V[end], Diagonal(ones(n))) # last covarianec matrix is identity

	# a tuple of m d-by-d variance component parameters
	Σ = ntuple(x -> zeros(d, d), m)
	for i in 1:m
	  copy!(Σ[i], generateSPDmatrix(d))
	end
	return(X, B, Σ, V)
end

import TraitSimulation: snparray_simulation
n = 10
p = 2
d = 2
m = 2

df = DataFrame(x = repeat([0.0], n), y = repeat([1.0], n))
dist = Normal()
link = IdentityLink()

# test for correct mean formula
formulas = ["x + 5y", "2 + log(y)"]

# what happens when there is no variables and just a scalar
formulas2 = ["25", "738"]
@test unique(mean_formula(formulas2[1], df)[1]) == [25]

evaluated_output = [repeat([5.0], n), repeat([2.0], n)]

for i in eachindex(formulas)
  return(@test mean_formula(formulas[i], df)[1] == evaluated_output[i])
end

X, B, Σ, V = generateRandomVCM(n, p, d, m)
test_vcm1 = VCMTrait(X, B, @vc Σ[1] ⊗ V[1] + Σ[2] ⊗ V[2])
test_vcm1_equivalent = VCMTrait(X, B, [Σ...], [V...])
@test test_vcm1_equivalent.vc[1].V == V[1]
@test typeof(test_vcm1.vc[1]) == VarianceComponent

varcomp = @vc Σ[1] ⊗ V[1] + Σ[2] ⊗ V[2]

varcomp_onevc = @vc Σ[1] ⊗ V[1]
@test eltype(varcomp_onevc) == VarianceComponent

test_vcm1 = VCMTrait(X, B, varcomp)

# check if the structure is correct
@test eltype(varcomp) == VarianceComponent

# check if returns the appropriate decomposition of the VarianceComponent type
@test vcobjtuple(varcomp)[1][1] == Σ[1]

# test provided simulate coefficients function
x = rand(n)
@test eltype(TraitSimulation.simulate_effect_size(x)) == Float64

effectsizes = rand(n)
our_names = ["sarah"; "janet"; "hua"; "eric"; "ken"; "jenny"; "ben"; "chris"; "juhyun"; "xinkai"]
whats_my_mean_formula = TraitSimulation.FixedEffectTerms(effectsizes, our_names)
data_frame_2 = DataFrame(ones(n, n), :auto)
rename!(data_frame_2, Symbol.(our_names))

@test unique(mean_formula(whats_my_mean_formula, data_frame_2)[1])[1] == sum(effectsizes)

@test_throws ErrorException TraitSimulation.__default_behavior(test_vcm1)

test_vcm1_new = VCMTrait(X, B, @vc Σ[1] ⊗ V[1] + Σ[2] ⊗ V[2])
test_vcm1_equiv_new = VCMTrait(X, B, [Σ...], [V...])

##
nsim = 10
using Statistics
Y_new = simulate(test_vcm1_new, nsim)
Y_vecd = zeros(n*d, nsim)

for i in 1:nsim
	Y_vecd[:, i]  = vec(Y_new[i])
end

simulated_mean = Statistics.mean(Y_vecd, dims = 2)

Z_new =  Y_vecd .- simulated_mean

emp_cov = (Z_new * Z_new') * inv(nsim)

true_mu = vec(test_vcm1_new.μ)

true_Ω = zeros(n*d, n*d)
for i = 1:m
  global true_Ω += kron(Σ[i], V[i])
end

vs = diag(true_Ω)

    for i = 1:20
        @test isapprox(simulated_mean[i], true_mu[i], atol=sqrt(vs[i] / nsim) * 8.0)
    end
    for i = 1:20, j = 1:20
        @test isapprox(emp_cov[i,j], true_Ω[i,j], atol=sqrt(vs[i] * vs[j]) * 10.0 / sqrt(nsim))
    end


text:
using Test, Hose

import Base: macroexpand
macroexpand(q) = macroexpand(Main, q)

#No change to nonhoses functionality
@test macroexpand( :(@hose a) ) == :a #doesn't change single inputs
@test macroexpand( :(@hose b(a)) ) == :(b(a)) #doesn't change inputs that a function applications

#Compatable with Julia 1.3 piping functionality
@test macroexpand( :(@hose a|>b) ) == :(b(a)) #basic
@test macroexpand( :(@hose a|>b|>c) ) == :(c(b(a)))  #Keeps chaining 3
@test macroexpand( :(@hose a|>b|>c|>d) ) == :(d(c(b(a)))) #Keeps chaining 4

@test macroexpand( :(@hose a |> b(x)) ) == :(b(a, x))  # Applying to function adds the LHS as the first argument Lazy.jl style
@test macroexpand( :(@hose a(x)|>b ) ) == :(b(a(x)))   #feeding functioncall results on wards

@test macroexpand(:(@hose 1|>a)) ==:(a(1)) #Works with literals (int)
@test macroexpand(:(@hose "foo"|>a)) == :(a("foo")) #Works with literal (string)
@test macroexpand( :(@hose a|>bb[2])) == :((bb[2])(a)) #Should work with RHS that is a array reference

# Test Problem
@test macroexpand( :(@hose a|>bb[2])) == :((bb[2])(a)) #Should work with RHS that is a array reference

#Marked locations
@test macroexpand( :(@hose a |> _)) == :(a) #Identity works
@test macroexpand( :(@hose a |> _[b])) == :(a[b]) #Indexing works

@test macroexpand( :(@hose a|>b(_) ) ) == :(b(a)) #Marked location only
@test macroexpand( :(@hose a|>b(x,_) ) ) == :(b(x,a)) # marked 2nd (and last)
@test macroexpand( :(@hose a|>b(_,x) ) ) == :(b(a,x)) # marked first
@test macroexpand( :(@hose a|>b(_,_) ) ) == :(b(a,a)) # marked double (Not certain if this is a good idea)
@test macroexpand( :(@hose a|>bb[2](x,_))) == :((bb[2])(x,a)) #Should work with RHS that is a array reference

#Macros and blocks
macro testmacro1(arg)
    esc(:($arg + 42))
end
macro testmacro2(arg, n)
    esc(:($arg + $n))
end
@test macroexpand( :(@hose a |> @testmacro1 ) ) == :(a + 42) # Can hose into macros like functions Lazy.jl style
@test macroexpand( :(@hose a |> @testmacro2 _ 3 ) ) == :(a + 3) # Can hose into macros
@test macroexpand( :(@hose a |> begin b = _; c + b + _ end )) == :(
                                begin b = a; c + b + a end)

#marked Unpacking
@test macroexpand( :(@hose a|>b(_...) ) ) == :(b(a...)) # Unpacking
@test macroexpand( :(@hose a|>bb[2](_...))) == :((bb[2])(a...)) #Should work with RHS of arry ref and do unpacking

#Mixing modes
@test macroexpand( :(@hose a|>b|>c(_) ) ) == :(c(b(a)))
@test macroexpand( :(@hose a|>b(x,_)|>c|>d(_,y) ) ) == :(d(c(b(x,a)),y))


text:
# example 2: read with readentry
#            default return type (String)

module FastaIOExample2
using FastaIO

function read_fasta_file(filename::String)
    # note: FastaReader can take a readable IO descriptor as
    #       argument instead of a filename; it also supports
    #       do-notation:
    #
    #       FastaReader(filename) do fr
    #         # do something
    #       end
    fr = FastaReader(filename)

    while !eof(fr)
        desc, seq = readentry(fr)
        println("num=$(fr.num_parsed) desc=$desc seq=$seq")
    end
    println("read $(fr.num_parsed) entries")
end

examples_dir = dirname(Base.source_path())
read_fasta_file(joinpath(examples_dir, "example.fasta.gz"))
end


text:
__precompile__(true)
module MathematicalSystems

using LinearAlgebra, SparseArrays
using LinearAlgebra: checksquare
using RecipesBase

#=======================
Identity operator
=======================#
include("identity.jl")

export IdentityMultiple, I

#=========================
Abstract Types for Systems
==========================#
include("abstract.jl")

#=======================
Utility functions
=======================#
include("utilities.jl")

# types
export AbstractSystem,
       AbstractDiscreteSystem,
       AbstractContinuousSystem

# methods
export statedim,
       inputdim,
       noisedim,
       stateset,
       inputset,
       noiseset,
       state_matrix,
       input_matrix,
       noise_matrix,
       affine_term,
       mass_matrix,
       viscosity_matrix,
       stiffness_matrix

# traits
export islinear,
       isaffine,
       ispolynomial,
       isblackbox,
       isnoisy,
       iscontrolled,
       isconstrained

#====================================
Concrete Types for Continuous Systems
====================================#
include("systems.jl")

export ContinuousIdentitySystem,
       ConstrainedContinuousIdentitySystem,
       LinearContinuousSystem,
       AffineContinuousSystem,
       LinearControlContinuousSystem,
       AffineControlContinuousSystem,
       ConstrainedLinearContinuousSystem,
       ConstrainedAffineContinuousSystem,
       ConstrainedAffineControlContinuousSystem,
       ConstrainedLinearControlContinuousSystem,
       LinearAlgebraicContinuousSystem,
       ConstrainedLinearAlgebraicContinuousSystem,
       PolynomialContinuousSystem,
       ConstrainedPolynomialContinuousSystem,
       BlackBoxContinuousSystem,
       ConstrainedBlackBoxContinuousSystem,
       BlackBoxControlContinuousSystem,
       ConstrainedBlackBoxControlContinuousSystem,
       NoisyLinearContinuousSystem,
       NoisyConstrainedLinearContinuousSystem,
       NoisyLinearControlContinuousSystem,
       NoisyConstrainedLinearControlContinuousSystem,
       NoisyAffineControlContinuousSystem,
       NoisyConstrainedAffineControlContinuousSystem,
       NoisyBlackBoxControlContinuousSystem,
       NoisyConstrainedBlackBoxControlContinuousSystem,
       SecondOrderLinearContinuousSystem,
       SecondOrderAffineContinuousSystem,
       SecondOrderConstrainedLinearControlContinuousSystem,
       SecondOrderConstrainedAffineControlContinuousSystem,
       SecondOrderContinuousSystem,
       SecondOrderConstrainedContinuousSystem

#==================================
Concrete Types for Discrete Systems
===================================#
export DiscreteIdentitySystem,
       ConstrainedDiscreteIdentitySystem,
       LinearDiscreteSystem,
       AffineDiscreteSystem,
       LinearControlDiscreteSystem,
       AffineControlDiscreteSystem,
       ConstrainedLinearDiscreteSystem,
       ConstrainedAffineDiscreteSystem,
       ConstrainedLinearControlDiscreteSystem,
       ConstrainedAffineControlDiscreteSystem,
       LinearAlgebraicDiscreteSystem,
       ConstrainedLinearAlgebraicDiscreteSystem,
       PolynomialDiscreteSystem,
       ConstrainedPolynomialDiscreteSystem,
       BlackBoxDiscreteSystem,
       ConstrainedBlackBoxDiscreteSystem,
       BlackBoxControlDiscreteSystem,
       ConstrainedBlackBoxControlDiscreteSystem,
       NoisyLinearDiscreteSystem,
       NoisyConstrainedLinearDiscreteSystem,
       NoisyLinearControlDiscreteSystem,
       NoisyConstrainedLinearControlDiscreteSystem,
       NoisyAffineControlDiscreteSystem,
       NoisyConstrainedAffineControlDiscreteSystem,
       NoisyBlackBoxControlDiscreteSystem,
       NoisyConstrainedBlackBoxControlDiscreteSystem,
       SecondOrderLinearDiscreteSystem,
       SecondOrderAffineDiscreteSystem,
       SecondOrderConstrainedLinearControlDiscreteSystem,
       SecondOrderConstrainedAffineControlDiscreteSystem,
       SecondOrderDiscreteSystem,
       SecondOrderConstrainedDiscreteSystem

#==========================================
Concrete Types for an Initial Value Problem
===========================================#
include("ivp.jl")

export InitialValueProblem, IVP,
       initial_state,
       system

#=====================
Input related methods
=====================#
include("inputs.jl")

export AbstractInput,
       ConstantInput,
       VaryingInput,
       nextinput

#==================================
Maps
===================================#
include("maps.jl")

# types
export AbstractMap,
       IdentityMap,
       ConstrainedIdentityMap,
       LinearMap,
       ConstrainedLinearMap,
       AffineMap,
       ConstrainedAffineMap,
       LinearControlMap,
       ConstrainedLinearControlMap,
       AffineControlMap,
       ConstrainedAffineControlMap,
       ResetMap,
       ConstrainedResetMap

# methods
export outputmap,
       outputdim,
       apply

#=========================
Systems with outputs
==========================#
include("outputs.jl")

export SystemWithOutput,
       LinearTimeInvariantSystem,
       LTISystem

#=========================
Macros
==========================#
include("macros.jl")

export @map,
       @system,
       @ivp

#===================================
Successor state for discrete systems
and vector field for continuous systems
====================================#
include("instantiate.jl")
include("successor.jl")
include("vector_field.jl")

export successor,
       vector_field

#===================================
Discretization for affine systems
====================================#
include("discretize.jl")

export discretize

end # module


text:
# include -
include("Include.jl")

# Load ensemble -
ensemble_array = readdlm("./simple_ensemble_parameter_sets.dat")

# What is the size of the ensemble?
[number_of_parameters,number_of_samples] = size(ensemble_array)

# Last ... run the model w/these sets, save the data -
for parameter_set_index = 1:number_of_samples

  # sample the ensemble_array -
  parameter_array = ensemble_array[:,parameter_set_index]

  # Load the data dictionary -
  data_dictionary = DataDictionary(0.0,0.0,0.0);

  # Update data dictionary to match new parameters before calculating obj
  parameter_mapping_array = data_dictionary["parameter_name_mapping_array"]
  for index = 1:length(parameter_mapping_array)
    if index <= data_dictionary["number_of_binding"]
      data_dictionary["binding_parameter_dictionary"][parameter_mapping_array[index]] = parameter_array[index]
    else
      data_dictionary["control_parameter_dictionary"][parameter_mapping_array[index]] = parameter_array[index]
    end
  end

  # Setup the time -
  time_start = 0.0
  time_stop = 120.0
  time_step_size = 0.01

  # Run the simulation -
  (T,X) = Simulation(time_start,time_stop,time_step_size,data_dictionary);

  # dump data to disk -
  local_data = [T X];
  data_filename = "./ensemble/sim_data_EI_"*string(parameter_set_index)*"_LI_"*string(local_index)*".dat"
  writedlm(data_filename,local_data);
end


text:
###  Some duas variáveis
# Faça um programa que some dois valores e exiba o resultado

print("inform a number ")
sum = parse(UInt8, readline())

print("inform other number ")
sum += parse(UInt8, readline())

print(sum,"\n")

max_stars_repo_name: codinginbrazil/functional-programming

max_issues_repo_name: codinginbrazil/functional-programming

max_forks_repo_name: codinginbrazil/functional-programming


text:
module Middle

include("./../../constraints/geometric/constants.jl")
include("./1_get_supporting.jl")
include("./2_incentive.jl")
using JuMP
using .Constants
using .GetSupporting
using .Incentive

export cons_ds_middle

function cons_ds_middle(m)

    m = mid_get_supporting(m)
    m = incentive(m)

    @expression(m, middle_objective,
        sum(m[:ds_middle_minmax_objective]) +
        sum(m[:ds_incentive_obj])
    )

    return m
end

end


text:
"""
    gumbelfitpwm(...)

Estimate the Gumbel parameters with the probability weighted moments.

# Implementation

Estimation with the probability weighted moments, as described by [Landwehr *et al. (1979)](https://agupubs.onlinelibrary.wiley.com/doi/abs/10.1029/WR015i005p01055),
is only possible in the stationary case.

See also [`gumbelfitpwm`](@ref) for the other methods, [`gevfit`](@ref), [`gevfitbayes`](@ref) and [`BlockMaxima`](@ref).

# Reference

Landwehr, J. M., Matalas, N. C. and Wallis, J. R. (1979). Probability weighted moments compared with
    some traditional techniques in estimating Gumbel parameters and quantiles. *Water Resources Research*,
    15:1055–1064.
"""
function gumbelfitpwm end

"""
    gumbelfitpwm(y::Vector{<:Real})

Estimate the Gumbel parameters with the probability weighted moments.

"""
function gumbelfitpwm(y::Vector{<:Real})::pwmEVA

    model = BlockMaxima(Variable("y", y))

    fittedmodel = gumbelfitpwm(model)

    return fittedmodel

end

"""
    gumbelfitpwm(df::DataFrame, datacol::Symbol)::pwmEVA

Estimate the Gumbel parameters with the probability weighted moments.

Block maxima data are in the column `datacol` of the dataframe `df`.
"""
function gumbelfitpwm(df::DataFrame, datacol::Symbol)::pwmEVA

    model = BlockMaxima(Variable(string(datacol), df[:, datacol]))

    fittedmodel = gumbelfitpwm(model)

    return fittedmodel

end

"""
    gumbelfitpwm(model::BlockMaxima)

Estimate the Gumbel parameters with the probability weighted moments.
"""
function gumbelfitpwm(model::BlockMaxima)::pwmEVA

    model = validatestationarity(model)

    y = model.data.value

    a₀ = pwm(y,1,0,0)
    a₁ = pwm(y,1,0,1)

    σ̂ = (a₀ - 2a₁)/log(2)

    μ̂ = a₀ - Base.MathConstants.eulergamma*σ̂

    ϕ̂ = log(σ̂)

    θ̂ = [μ̂, ϕ̂, 0.0]

    fm = pwmEVA{BlockMaxima, Gumbel}(model,θ̂)

    return fm

end


text:
## RESAMPLING STRATEGIES

abstract type ResamplingStrategy <: MLJType end
show_as_constructed(::Type{<:ResamplingStrategy}) = true

# resampling strategies are `==` if they have the same type and their
# field values are `==`:
function ==(s1::S, s2::S) where S <: ResamplingStrategy
    return all(getfield(s1, fld) == getfield(s2, fld) for fld in fieldnames(S))
end

train_test_pairs(s::ResamplingStrategy, rows, X, y, w) =
    train_test_pairs(s, rows, X, y)
train_test_pairs(s::ResamplingStrategy, rows, X, y) =
    train_test_pairs(s, rows)

"""
    Holdout(; fraction_train=0.7,
              shuffle=false,
              rng=Random.GLOBAL_RNG)

Single train-test split with a (randomly selected) portion of the
data being selected for training and the rest for testing.

If `rng` is an integer, then `MersenneTwister(rng)` is the random
number generator used for shuffling rows. Otherwise some `AbstractRNG`
object is expected.

"""
struct Holdout <: ResamplingStrategy
    fraction_train::Float64
    shuffle::Bool
    rng::Union{Int,AbstractRNG}

    function Holdout(fraction_train, shuffle, rng)
        0 < fraction_train < 1 ||
            error("`fraction_train` must be between 0 and 1.")
        return new(fraction_train, shuffle, rng)
    end
end

# Keyword Constructor
function Holdout(; fraction_train::Float64=0.7,
                   shuffle::Bool=false,
                   rng::Union{Int,AbstractRNG}=Random.GLOBAL_RNG)
    Holdout(fraction_train, shuffle, rng)
end

function train_test_pairs(holdout::Holdout, rows)
    if holdout.rng isa Integer
        rng = MersenneTwister(holdout.rng)
    else
        rng = holdout.rng
    end
    train, test = partition(rows, holdout.fraction_train,
                             shuffle=holdout.shuffle, rng=rng)
    return [(train, test),]
end

"""
    CV(; nfolds=6,  shuffle=false, rng=Random.GLOBAL_RNG)

Cross validation resampling where the data is (randomly) partitioned
in `nfolds` folds and the model is evaluated `nfolds` times, each time
taking one fold for testing and the remaining folds for training.

For instance, if `nfolds=3` then the data will be partitioned in three
folds A, B and C and the model will be trained three times, first with
A and B and tested on C, then on A, C and tested on B and finally on
B, C and tested on A.

If `rng` is an integer, then `MersenneTwister(rng)` is the random
number generator used for shuffling rows. Otherwise some `AbstractRNG`
object is expected.

"""
struct CV <: ResamplingStrategy
    nfolds::Int
    shuffle::Bool
    rng::Union{Int,AbstractRNG}
    function CV(nfolds, shuffle, rng)
        nfolds > 1 || error("Must have nfolds > 1. ")
        return new(nfolds, shuffle, rng)
    end
end

# Constructor with keywords
CV(; nfolds::Int=6,  shuffle::Bool=false,
   rng::Union{Int,AbstractRNG}=Random.GLOBAL_RNG) =
       CV(nfolds, shuffle, rng)

function train_test_pairs(cv::CV, rows)
    if cv.rng isa Integer
        rng = MersenneTwister(cv.rng)
    else
        rng = cv.rng
    end

    n_observations = length(rows)
    nfolds = cv.nfolds

    if cv.shuffle
        rows=shuffle!(rng, collect(rows))
    end

    # number of observations per fold
    k = floor(Int, n_observations/nfolds)
    k > 0 || error("Inusufficient data for $nfolds-fold cross-validation.\n"*
                   "Try reducing nfolds. ")

    # define the (trainrows, testrows) pairs:
    firsts = 1:k:((nfolds - 1)*k + 1) # itr of first `test` rows index
    seconds = k:k:(nfolds*k)          # itr of last  `test` rows index
    ret = map(1:nfolds) do k
        f = firsts[k]
        s = seconds[k]
        return (vcat(rows[1:(f - 1)], rows[(s + 1):end]), # trainrows
                rows[f:s])                               # testrows
    end

    return ret
end

## DIRECT EVALUATION METHODS

"""
    evaluate!(mach,
              resampling=CV(),
              measure=nothing,
              weights=nothing,
              operation=predict,
              acceleration=DEFAULT_RESOURCE[],
              force=false,
              verbosity=1)

Estimate the performance of a machine `mach` wrapping a supervised
model in data, using the specified `resampling` strategy (defaulting
to 6-fold cross-validation) and `measure`, which can be a single
measure or vector.

Do `subtypes(MLJ.ResamplingStrategy)` to obtain a list of available resampling
strategies. If `resampling` is not an object of type
`MLJ.ResamplingStrategy`, then a vector of pairs (of the form
`(train_rows, test_rows)` is expected. For example, setting

    resampling = [(1:100), (101:200)),
                   (101:200), (1:100)]

gives two-fold cross-validation using the first 200 rows of data.

If `resampling isa MLJ.ResamplingStrategy` then one may optionally
restrict the data used in evaluation by specifying `rows`.

An optional `weights` vector may be passed for measures that support
sample weights (`MLJ.supports_weights(measure) == true`), which
is ignored by those that don't.

User-defined measures are supported; see the manual for details.

If no measure is specified, then `default_measure(mach.model)` is
used, unless this default is `nothing` and an error is thrown.

The `acceleration` keyword argument is used to specify the compute resource (a
subtype of `ComputationalResources.AbstractResource`) that will be used to
accelerate/parallelize the resampling operation.

Although evaluate! is mutating, `mach.model` and `mach.args` are
untouched.

"""
evaluate!(mach::Machine{<:Supervised};
          resampling=CV(), kwargs...) =
              evaluate!(mach, resampling; kwargs... )

"""
    evaluate(model, X, y; measure=nothing, options...)

Evaluate the performance of a supervised model `model` on input data
`X` and target `y`. See the machine version `evaluate!` for options.

"""
evaluate(model::Supervised, args...; kwargs...) =
    evaluate!(machine(model, args...); kwargs...)

const AbstractRow = Union{AbstractVector{<:Integer}, Colon}
const TrainTestPair = Tuple{AbstractRow,AbstractRow}
const TrainTestPairs = AbstractVector{<:TrainTestPair}

function _check_measure(model, measure, y, operation, override)

    override && (return nothing)

    T = scitype(y)

    T == Unknown && (return nothing)
    target_scitype(measure) == Unknown && (return nothing)
    prediction_type(measure) == :unknown && (return nothing)

    avoid = "\nTo override measure checks, set check_measure=false. "

    T <: target_scitype(measure) ||
        throw(ArgumentError(
            "\nscitype of target = $T but target_scitype($measure) = "*
            "$(target_scitype(measure))."*avoid))

    if model isa Probabilistic
        if operation == predict
            if prediction_type(measure) != :probabilistic
                suggestion = ""
                if target_scitype(measure) <: Finite
                    suggestion = "\nPerhaps you want to set operation="*
                    "predict_mode. "
                elseif target_scitype(measure) <: Continuous
                    suggestion = "\nPerhaps you want to set operation="*
                    "predict_mean or operation=predict_median. "
                else
                    suggestion = ""
                end
                throw(ArgumentError(
                    "\n$model <: Probabilistic but prediction_type($measure) = "*
                      ":$(prediction_type(measure)). "*suggestion*avoid))
            end
        end
    end

    model isa Deterministic && prediction_type(measure) != :deterministic &&
        throw(ArgumentError("$model <: Deterministic but prediction_type($measure) ="*
              ":$(prediction_type(measure))."*avoid))

    return nothing

end

function _evaluate!(func::Function, res::CPU1, nfolds, verbosity)
    p = Progress(nfolds + 1, dt=0, desc="Evaluating over $nfolds folds: ",
                 barglyphs=BarGlyphs("[=> ]"), barlen=25, color=:yellow)
    verbosity > 0 && next!(p)
    return reduce(vcat, (func(k, p, verbosity) for k in 1:nfolds))
end
function _evaluate!(func::Function, res::CPUProcesses, nfolds, verbosity)
    # TODO: use pmap here ?:
    return @distributed vcat for k in 1:nfolds
        func(k)
    end
end
@static if VERSION >= v"1.3.0-DEV.573"
    function _evaluate!(func::Function, res::CPUThreads, nfolds, verbosity)
        task_vec = [Threads.@spawn func(k) for k in 1:nfolds]
        return fetch.(task_vec)
    end
end

# if `resampling` is not a ResamplingStrategy object:
function evaluate!(mach::Machine, resampling;
                   measure=nothing, weights=nothing,
                   operation=predict, acceleration=DEFAULT_RESOURCE[],
                   rows=nothing, force=false,
                   check_measure=true, verbosity=1)

    resampling isa TrainTestPairs ||
        error("`resampling` must be an "*
              "MLJ.ResamplingStrategy or tuple of pairs "*
              "of the form `(train_rows, test_rows)`")

    rows === nothing ||
        error("You cannot specify `rows` unless `resampling "*
              "isa MLJ.ResamplingStrategy`. ")

    if measure === nothing
        candidate = default_measure(mach.model)
        candidate ===  nothing && error("You need to specify measure=... ")
        measures = [candidate, ]
    elseif !(measure isa AbstractVector)
        measures = [measure, ]
    else
        measures = measure
    end

    X = mach.args[1]
    y = mach.args[2]

    [_check_measure(mach.model, m, y, operation, !check_measure) for m in measures]

    if verbosity >= 0 && weights !== nothing
        unsupported = filter(measures) do m
            !supports_weights(m)
        end
        if !isempty(unsupported)
            unsupported_as_string = string(unsupported[1])
            unsupported_as_string *=
                reduce(*, [string(", ", m) for m in unsupported[2:end]])
            @warn "weights ignored in evaluations of the following measures, "*
            "as unsupported: \n$unsupported_as_string "
        end
    end

    nfolds = length(resampling)
    nmeasures = length(measures)

    function get_measurements(k)
        train, test = resampling[k]
        fit!(mach; rows=train, verbosity=verbosity-1, force=force)
        Xtest = selectrows(X, test)
        ytest = selectrows(y, test)
        if weights == nothing
            wtest = nothing
        else
            wtest = weights[test]
        end
        yhat = operation(mach, Xtest)
        return [value(m, yhat, Xtest, ytest, wtest)
                for m in measures]
    end
    function get_measurements(k, p, verbosity) # p = progress meter
        ret = get_measurements(k)
        verbosity > 0 && next!(p)
        return ret
    end

    measurements_flat = if acceleration isa CPUProcesses
        ## TODO: progress meter for distributed case
        if verbosity > 0
            @info "Distributing cross-validation computation " *
                  "among $(nworkers()) workers."
        end
    end
    measurements_flat = _evaluate!(get_measurements, acceleration, nfolds, verbosity)

    # in the following rows=folds, columns=measures:
    measurements_matrix = permutedims(
        reshape(measurements_flat, (nmeasures, nfolds)))

    # measurements for each observation:
    per_observation = map(1:nmeasures) do k
        m = measures[k]
        if reports_each_observation(m)
            [measurements_matrix[:,k]...]
        else
            missing
        end
    end

    # measurements for each fold:
    per_fold = map(1:nmeasures) do k
        m = measures[k]
        if reports_each_observation(m)
            broadcast(MLJBase.aggregate, per_observation[k], [m,])
        else
            [measurements_matrix[:,k]...]
        end
    end

    # overall aggregates:
    per_measure = map(1:nmeasures) do k
        m = measures[k]
        MLJBase.aggregate(per_fold[k], m)
    end

    ret = (measure=measures,
           measurement=per_measure,
           per_fold=per_fold,
           per_observation=per_observation)

    verbosity < 1 || nmeasures < 2 ||
        pretty(selectcols(ret, 1:2), showtypes=false)

    return ret

end

function actual_rows_and_weights(rows, weights, N, verbosity)
    unspecified_rows = (rows === nothing)
    _rows = unspecified_rows ? (1:N) : rows
    if weights == nothing
        _weights = nothing
    else
        _weights = weights[_rows]
    end
    unspecified_rows || @info "Creating subsamples from a subset of all rows. "
    return _rows, _weights
end

# evaluation when ResamplingStrategy is passed (instead of train/test rows):
function evaluate!(mach::Machine, resampling::ResamplingStrategy;
                   weights=nothing, rows=nothing, verbosity=1, kwargs...)

    y = mach.args[2]
    _rows, _weights =
        actual_rows_and_weights(rows, weights, length(y), verbosity)

    return evaluate!(mach::Machine,
               train_test_pairs(resampling, _rows, mach.args...);
               weights=_weights, verbosity=verbosity, kwargs...)

end

## RESAMPLER - A MODEL WRAPPER WITH `evaluate` OPERATION

"""
$TYPEDEF

Resampler structure for the `TunedModel` `fit` defined in `tuning.jl`.
"""
mutable struct Resampler{S,M<:Supervised} <: Supervised
    model::M
    resampling::S # resampling strategy
    measure
    weights::Union{Nothing,AbstractVector{<:Real}}
    operation
    acceleration::AbstractResource
end

MLJBase.package_name(::Type{<:Resampler}) = "MLJ"
MLJBase.is_wrapper(::Type{<:Resampler}) = true

Resampler(; model=ConstantRegressor(), resampling=Holdout(),
            measure=nothing, weights=nothing, operation=predict,
            acceleration=DEFAULT_RESOURCE[]) =
                Resampler(model, resampling, measure, weights, operation,
                          acceleration)

function MLJBase.fit(resampler::Resampler, verbosity::Int, X, y)

    if resampler.measure === nothing
        measure = default_measure(resampler.model)
        if measure === nothing
            error("You need to specify measure=... ")
        end
    else
        measure = resampler.measure
    end

    mach = machine(resampler.model, X, y)

    fitresult = evaluate!(mach, resampler.resampling;
                          measure=measure, weights=resampler.weights,
                          operation=resampler.operation,
                          verbosity=verbosity-1,
                          acceleration=resampler.acceleration)

    cache = (mach, deepcopy(resampler.resampling))
    report = NamedTuple()

    return fitresult, cache, report

end

# in special case of holdout, we can reuse the underlying model's
# machine, provided the training_fraction has not changed:
function MLJBase.update(resampler::Resampler{Holdout},
                        verbosity::Int, fitresult, cache, X, y)

    old_mach, old_resampling = cache

    if resampler.measure === nothing
        measure = default_measure(resampler.model)
        if measure === nothing
            error("You need to specify measure=... ")
        end
    else
        measure = resampler.measure
    end

    if old_resampling.fraction_train == resampler.resampling.fraction_train
        mach = old_mach
    else
        mach = machine(resampler.model, X, y)
        cache = (mach, deepcopy(resampler.resampling))
    end

    fitresult = evaluate!(mach, resampler.resampling;
                          measure=resampler.measure,
                          weights=resampler.weights,
                          operation=resampler.operation,
                          verbosity=verbosity-1)

    report = NamedTuple

    return fitresult, cache, report

end

MLJBase.input_scitype(::Type{<:Resampler{S,M}}) where {S,M} = MLJBase.input_scitype(M)
MLJBase.target_scitype(::Type{<:Resampler{S,M}}) where {S,M} = MLJBase.target_scitype(M)

## EVALUATE FOR MODEL + DATA

MLJBase.evaluate(model::Resampler, fitresult) = fitresult


text:
# function tcm_estimateForest_kcv(N, d_list, nfolds, transaction_folds, transactions, time_limit, perTol, initial_left, initial_right, initial_product, initial_isLeaf)
# 	unique_fold_ids = unique(transaction_folds)

# 	KL_avg_holdout_fold = zeros(length(d_list), nfolds)

# 	start_time = time()

# 	for fold_ind in 1:nfolds
# 		fold = unique_fold_ids[fold_ind]

# 		@show fold

# 		@show unique_fold_ids

# 		println("Pre-subsetting ")

# 		#@show transaction_folds
# 		#temp = find(transaction_folds .!= fold)
# 		#@show temp
# 		@show size(transaction_folds)
# 		@show size(transactions)
# 		@show fold

# 		temp = find(transaction_folds .!= fold)
# 		@show temp[1:5]
# 		trainfolds_transactions = transactions[ find(transaction_folds .!= fold), :]
# 		println("Post-subsetting")
# 		testfold_transactions = transactions[ find(transaction_folds .== fold), :]

# 		trainfolds_assortments, trainfolds_transaction_counts = tcm_transactionsToCounts(N, trainfolds_transactions)
# 		testfold_assortments, testfold_transaction_counts = tcm_transactionsToCounts(N, testfold_transactions)

# 		for i in 1:length(d_list)
# 			println("Depth ", d_list[i], ", fold ", fold)
# 			depth_limit = d_list[i]

# 			lambda, forest_left, forest_right, forest_product, forest_isLeaf, loglik, elapsed_time = tcm_estimateForest(N, trainfolds_assortments, trainfolds_transaction_counts, time_limit, perTol, depth_limit,
# 																												initial_left, initial_right, initial_product, initial_isLeaf)
# 			K = length(lambda)
# 			testfold_forest_predict, testfold_forest_predict_by_class = tcm_predictForest(N, K, lambda, forest_left, forest_right, forest_product, forest_isLeaf, testfold_assortments)

# 			testfold_KL, testfold_noinf_count = tcm_evaluateKL(N, testfold_forest_predict, testfold_transaction_counts)

# 			KL_avg_holdout_fold[i, fold_ind] = testfold_KL / testfold_noinf_count #testfold_KL / sum(testfold_transaction_counts)
# 		end
# 	end

# 	end_time = time()
# 	elapsed_time = end_time - start_time

# 	KL_avg_by_d = mean(KL_avg_holdout_fold, 2) * 100

# 	best_i = indmin(KL_avg_by_d)
# 	best_d = d_list[best_i]
# 	best_KL_avg_cv = KL_avg_by_d[best_i]

# 	return best_d, best_KL_avg_cv, elapsed_time

# end

function tcm_estimateForest_kcv_weighted(N, d_list, nfolds, transaction_folds, transactions, time_limit, perTol, initial_left, initial_right, initial_product, initial_isLeaf)
	unique_fold_ids = unique(transaction_folds)

	KL_total_holdout_fold = zeros(length(d_list), nfolds)
	noinf_count_holdout_fold = zeros(length(d_list), nfolds)

	start_time = time()

	for fold_ind in 1:nfolds
		fold = unique_fold_ids[fold_ind]

		@show fold

		@show unique_fold_ids

		println("Pre-subsetting ")

		#@show transaction_folds
		#temp = find(transaction_folds .!= fold)
		#@show temp
		@show size(transaction_folds)
		@show size(transactions)
		@show fold

		temp = find(transaction_folds .!= fold)
		@show temp[1:5]
		trainfolds_transactions = transactions[ find(transaction_folds .!= fold), :]
		println("Post-subsetting")
		testfold_transactions = transactions[ find(transaction_folds .== fold), :]

		trainfolds_assortments, trainfolds_transaction_counts = tcm_transactionsToCounts(N, trainfolds_transactions)
		testfold_assortments, testfold_transaction_counts = tcm_transactionsToCounts(N, testfold_transactions)

		for i in 1:length(d_list)
			println("Depth ", d_list[i], ", fold ", fold)
			depth_limit = d_list[i]

			lambda, forest_left, forest_right, forest_product, forest_isLeaf, loglik, elapsed_time = tcm_estimateForest(N, trainfolds_assortments, trainfolds_transaction_counts, time_limit, perTol, depth_limit,
																												initial_left, initial_right, initial_product, initial_isLeaf)
			K = length(lambda)
			testfold_forest_predict, testfold_forest_predict_by_class = tcm_predictForest(N, K, lambda, forest_left, forest_right, forest_product, forest_isLeaf, testfold_assortments)

			testfold_KL, testfold_noinf_count, testfold_KL_avg = tcm_evaluateKL(N, testfold_forest_predict, testfold_transaction_counts)

			KL_total_holdout_fold[i, fold_ind] = testfold_KL #testfold_KL / sum(testfold_transaction_counts)
			noinf_count_holdout_fold[i, fold_ind] = testfold_noinf_count
		end
	end

	end_time = time()
	elapsed_time = end_time - start_time

	KL_avg_by_d = sum(KL_total_holdout_fold, 2) ./ sum(noinf_count_holdout_fold, 2)  * 100

	@show KL_total_holdout_fold
	@show noinf_count_holdout_fold
	@show KL_avg_by_d

	best_i = indmin(KL_avg_by_d)
	best_d = d_list[best_i]
	best_KL_avg_cv = KL_avg_by_d[best_i]

	return best_d, best_KL_avg_cv, elapsed_time

end

function tcm_estimateForest_kcv_weighted_ass(N, d_list, nfolds, assortment_folds, assortments, transaction_counts, time_limit, perTol, initial_left, initial_right, initial_product, initial_isLeaf)
	unique_fold_ids = unique(assortment_folds)

	KL_total_holdout_fold = zeros(length(d_list), nfolds)
	noinf_count_holdout_fold = zeros(length(d_list), nfolds)

	start_time = time()

	for fold_ind in 1:nfolds
		fold = unique_fold_ids[fold_ind]

		# trainfolds_assortments, trainfolds_transaction_counts = tcm_transactionsToCounts(N, trainfolds_transactions)
		# testfold_assortments, testfold_transaction_counts = tcm_transactionsToCounts(N, testfold_transactions)

		trainfolds_assortments = assortments[ find(assortment_folds .!= fold), :]
		testfold_assortments = assortments[ find(assortment_folds .== fold), :]

		trainfolds_transaction_counts = transaction_counts[ find(assortment_folds .!= fold), :]
		testfold_transaction_counts = transaction_counts[ find(assortment_folds .== fold), :]

		for i in 1:length(d_list)
			println("Depth ", d_list[i], ", fold ", fold)
			depth_limit = d_list[i]

			lambda, forest_left, forest_right, forest_product, forest_isLeaf, loglik, elapsed_time = tcm_estimateForest(N, trainfolds_assortments, trainfolds_transaction_counts, time_limit, perTol, depth_limit,
																												initial_left, initial_right, initial_product, initial_isLeaf)
			K = length(lambda)
			testfold_forest_predict, testfold_forest_predict_by_class = tcm_predictForest(N, K, lambda, forest_left, forest_right, forest_product, forest_isLeaf, testfold_assortments)

			testfold_KL, testfold_noinf_count, testfold_KL_avg = tcm_evaluateKL(N, testfold_forest_predict, testfold_transaction_counts)

			KL_total_holdout_fold[i, fold_ind] = testfold_KL #testfold_KL / sum(testfold_transaction_counts)
			noinf_count_holdout_fold[i, fold_ind] = testfold_noinf_count
		end
	end

	end_time = time()
	elapsed_time = end_time - start_time

	KL_avg_by_d = sum(KL_total_holdout_fold, 2) ./ sum(noinf_count_holdout_fold, 2)  * 100

	best_i = indmin(KL_avg_by_d)
	best_d = d_list[best_i]
	best_KL_avg_cv = KL_avg_by_d[best_i]

	return best_d, best_KL_avg_cv, elapsed_time

end

function tcm_estimateForest_RMredo_kcv_weighted_ass(N, d_list, nfolds, assortment_folds, assortments, transaction_counts, time_limit, perTol)
	unique_fold_ids = unique(assortment_folds)

	KL_total_holdout_fold = zeros(length(d_list), nfolds)
	noinf_count_holdout_fold = zeros(length(d_list), nfolds)

	start_time = time()

	for fold_ind in 1:nfolds
		fold = unique_fold_ids[fold_ind]

		# trainfolds_assortments, trainfolds_transaction_counts = tcm_transactionsToCounts(N, trainfolds_transactions)
		# testfold_assortments, testfold_transaction_counts = tcm_transactionsToCounts(N, testfold_transactions)

		trainfolds_assortments = assortments[ find(assortment_folds .!= fold), :]
		testfold_assortments = assortments[ find(assortment_folds .== fold), :]

		trainfolds_transaction_counts = transaction_counts[ find(assortment_folds .!= fold), :]
		testfold_transaction_counts = transaction_counts[ find(assortment_folds .== fold), :]

		for i in 1:length(d_list)
			println("Depth ", d_list[i], ", fold ", fold)
			depth_limit = d_list[i]

			ranking_lambda, orderings, ranking_loglik, ranking_elapsed_time = tcm_estimateRanking(N, trainfolds_assortments, trainfolds_transaction_counts, time_limit, perTol)
			initial_left, initial_right, initial_product, initial_isLeaf = tcm_convertRankingToForest(N, orderings)

			lambda, forest_left, forest_right, forest_product, forest_isLeaf, loglik, elapsed_time = tcm_estimateForest(N, trainfolds_assortments, trainfolds_transaction_counts, time_limit, perTol, depth_limit,
																												initial_left, initial_right, initial_product, initial_isLeaf)
			K = length(lambda)
			testfold_forest_predict, testfold_forest_predict_by_class = tcm_predictForest(N, K, lambda, forest_left, forest_right, forest_product, forest_isLeaf, testfold_assortments)

			testfold_KL, testfold_noinf_count, testfold_KL_avg = tcm_evaluateKL(N, testfold_forest_predict, testfold_transaction_counts)

			KL_total_holdout_fold[i, fold_ind] = testfold_KL #testfold_KL / sum(testfold_transaction_counts)
			noinf_count_holdout_fold[i, fold_ind] = testfold_noinf_count
		end
	end

	end_time = time()
	elapsed_time = end_time - start_time

	KL_avg_by_d = sum(KL_total_holdout_fold, 2) ./ sum(noinf_count_holdout_fold, 2)  * 100

	best_i = indmin(KL_avg_by_d)
	best_d = d_list[best_i]
	best_KL_avg_cv = KL_avg_by_d[best_i]

	return best_d, best_KL_avg_cv, elapsed_time

end

max_stars_repo_licenses: MIT Unlicense

max_issues_repo_licenses: MIT Unlicense

max_forks_repo_licenses: MIT Unlicense


text:
# Implementation of Karemakers adaption of the DeBoer Model
"""
    KaremakerModel()

A simple structure that stores a Karemaker-model.

The default initialization was taken from the original paper. When constructiong the model, single properties can be defined as keyword arguments.

# Examples
```julia
julia> model = KaremakerModel()
KaremakerModel
  c: Float64 350.0
  magn_int: Float64 10.0
  r_int: Distributions.Normal{Float64}
  A_resp: Float64 60.0
  T_resp: Float64 3000.0
  a1: Float64 9.0
  a2: Float64 9.0
  b_I: Float64 0.125
  b_D: Float64 0.125
  b_P: Float64 0.125
  Dref: Float64 80.0
  t_runoff: Float64 2000.0
  offset: Float64 10.0
  starling: Float64 0.016
  magn_pulse: Float64 2.0
  A_BP: Float64 3.0
  w: Array{Float64}((6,)) [0.0, 0.01, 0.03, 0.03, 0.015, 0.005]
  symp_set: Float64 0.5
  symp_add: Array{Float64}((2,)) [-0.05, 0.2]
```
"""
@with_kw mutable struct KaremakerModel
    c::Float64 = 350 # ms
    magn_int::Float64 = 10 # ms
    # r_int::Distribution = TriangularDist(-1, 1)
    r_int::Distribution = Normal(0, 1)
    # Respiration
    A_resp::Float64 = 60 # ms
    T_resp::Float64 = 3000 # ms
    # Baro
    a1::Float64 = 9 # mmHg/ms
    a2::Float64 = 9 # mmHg/ms
    # sympatic force in I
    b_I::Float64 = 0.125
    # sympatic force in D
    b_D::Float64 = 0.125
    # sympatic force in P
    b_P::Float64 = 0.125
    Dref::Float64 = 80 # mmHg
    t_runoff::Float64 = 2000 # ms
    offset::Float64 = 10 # mmHg
    starling::Float64 = 0.016 # mmHg/ms
    magn_pulse::Float64 = 2 # ms
    A_BP::Float64 = 3 # mmHg
    w::Vector{Float64} = [0.000,0.010, 0.030, 0.030, 0.015, 0.005]
    # w::Vector{Float64} = [1,2,3,3,2,1]
    # w::Vector{Float64} = [1,1,1,1,1,1]+
    ## Experimental ##
    symp_set = 0.5
    # symp_add = [-0.0, 0.0]
    symp_add = [-0.05, 0.2]
end

"""
    predict(model::KaremakerModel, n::Int; burnIn::Int = 0)

Predicts 'N' values of the cardiovascular variables for the respective model. With 'burnIn' a certain number af values can be dropped in the beginning.

# Examples
```julia
 julia> S, D, P, I, Symp, ρ = predict(model, 100)
 (S = [...], D = [...], P = [...], I = [...], Symp = [...], ρ = [...])
 ```
 """
function predict(model::KaremakerModel, n::Int; burnIn::Int = 50)
    t = Vector{Float64}(undef, n + burnIn)
    trm = Vector{Float64}(undef, n + burnIn)
    I = Vector{Float64}(undef, n + burnIn)
    D = Vector{Float64}(undef, n + burnIn)
    S = Vector{Float64}(undef, n + burnIn)
    P = Vector{Float64}(undef, n + burnIn)
    Symp = Vector{Float64}(undef, n + burnIn)
    SympD = Vector{Float64}(undef, n + burnIn)

    history = 6
    for i in 1:history
        t[i] = 0
        trm[i] = model.c + model.magn_int * rand(model.r_int) + model.A_resp * sin((2π * t[i]) / model.T_resp + 1.5π)
        I[i] = 1000
        D[i] = model.Dref
        S[i] = 160
        P[i] = 40
        Symp[i] = 0
        SympD[i] = 0

    end

    for i in (history+1):(n+burnIn)
        #timebase
        t[i] = t[i-1] + I[i-1]
        # constant + noise + rspiration for I
        trm[i] = model.c + model.magn_int * rand(model.r_int) + model.A_resp * sin((2π * t[i]) / model.T_resp + 1.5π)

        # Symp[i] = 2/π * atan((sum(model.w .* SympD[i .- (1:length(model.w))]))/ (sum(model.w)))
        # if Symp[i] > model.symp_set
        #     Symp[i] += model.symp_add[1] * Symp[i]
        # else
        #     Symp[i] += model.symp_add[2] * Symp[i]
        # end

        # normalize weights or not?
        Symp[i] = 2/π * atan((sum(model.w .* SympD[i .- (1:length(model.w))]))/ sum(model.w))
        # Symp[i] = atan((sum(model.w .* SympD[i .- (1:length(model.w))]))) # * 2/π
        if I[i-1] > 700
            I[i] = ((model.a1 * S[i-1] + trm[i]) / (1 + model.b_I * Symp[i]))
        else
            I[i] = ((model.a2 * S[i-2] + trm[i]) / (1 + model.b_I * Symp[i]))
        end
        D[i] = S[i-1] * exp(-I[i] / (model.t_runoff * (1 + model.b_D * Symp[i])))
        P[i] = (model.offset + model.starling * I[i]) * (1 + model.b_P * Symp[i]) + model.magn_pulse * rand(model.r_int) + model.A_BP * sin((2π * t[i]) / model.T_resp)
        S[i] = D[i] + P[i]
        if D[i] > model.Dref
            SympD[i] = 0
        else
            SympD[i] = model.Dref - D[i]
        end
    end
    ρ = model.A_resp .* sin.((2π .* t) ./ model.T_resp .+ 1.5π)
    return ( S = S[burnIn+1:end], D = D[burnIn+1:end], P = P[burnIn+1:end], I = I[burnIn+1:end], Symp = Symp[burnIn+1:end], ρ = ρ[burnIn+1:end])# , SympD = SympD[burnIn+1:end])
    # return (I = I[burnIn+1:end], S = S[burnIn+1:end], D = D[burnIn+1:end], Symp = Symp[burnIn+1:end])
end

# model = KaremakerModel()
# i, s, d, symp = predict(model, 1000, burnIn = 45)
# get_lfhf(i)
# mean(symp)
# pi = plot(i)
# plot!(pi, i)
# ps = plot(s)
# plot!(ps, s)
# pd = plot(d)
# plot!(pd, d)
# lfhfs = []
# symps = []
# rr = []
# for runoff in 1250:2850
#     model = KaremakerModel(t_runoff = runoff)
#     i, s, d, symp = predict(model, 1000, burnIn = 50)
#     pgram = welch_pgram(i)
#     inds_lf = findall(x -> 0.04 <= x <= 0.15, pgram.freq)
#     inds_hf = findall(x -> 0.15 <= x <= 0.4, pgram.freq)
#     pow_lf = sum(pgram.power[inds_lf])
#     pow_hf = sum(pgram.power[inds_hf])
#     lfhf = pow_lf / pow_hf
#     push!(lfhfs, lfhf)
#     push!(lfhfs, lfhf)
#     push!(symps, mean(symp))
#     push!(rr, mean(i))
# end

# plot(symps, rr)

function get_lfhf(i)
    fs = 1000 / (mean(i))
    i = i .- mean(i)
    pgram = welch_pgram(i, fs = fs)
    inds_lf = findall(x -> 0.04 <= x <= 0.15, pgram.freq)
    inds_hf = findall(x -> 0.15 <= x <= 0.5, pgram.freq)
    pow_lf = sum(pgram.power[inds_lf])
    pow_hf = sum(pgram.power[inds_hf])
    lfhf = pow_lf / pow_hf
    return lfhf
end

"""
    function vagal_balance(Dref::Real, runoffs::Union{AbstractVector{Real}, StepRangeLen}; n::Int = 5000, model::Union{Nothing, KaremakerModel} = nothing)

Calculates the LF/HF ratio for all distolic runoff constants in 'runoffs' returns also some means for different signals.
"""
function vagal_balance(Dref::Real, runoffs::Union{AbstractVector{Real}, StepRangeLen}; n::Int = 5000, model::Union{Nothing, KaremakerModel} = nothing)
    rr = []
    symp_drive = []
    lfhf = []
    sb = []
    db = []
    for runoff in runoffs
        if model  === nothing
            model = KaremakerModel(Dref = Dref, t_runoff = runoff)
        else
            model.Dref = Dref
            model.t_runoff = runoff
        end
        s, d, _, i, symp, _ = predict(model, n)
        # i, s, d, symp = predict(model, n)
        push!(rr, mean(i))
        push!(sb, mean(s))
        push!(db, mean(d))
        push!(lfhf, get_lfhf(i))
        push!(symp_drive, mean(symp))
    end
    return rr, symp_drive, lfhf, sb, db
end

function vagal_balance_respiration(Ts::Union{AbstractVector{Real}, StepRangeLen}, Dref::Real, runoffs::Union{AbstractVector{Real}, StepRangeLen}; n::Int = 5000, model::Union{Nothing, KaremakerModel} = nothing)
    rr = []
    symp_drive = []
    lfhf = []
    for T in Ts
        if model  === nothing
            model = KaremakerModel(T_resp = T)
        else
            model.T_resp = T
        end
        rr_T, symp_drive_T, lfhf_T = vagal_balance(Dref, runoffs, model = model)
        push!(rr, rr_T)
        push!(lfhf, lfhf_T)
        push!(symp_drive, symp_drive_T)
    end
    return rr, symp_drive, lfhf
end

function lfhf2D_respiration(Ts::Union{AbstractVector{Real}, StepRangeLen}, Dref::Real, runoffs::Union{AbstractVector{Real}, StepRangeLen}; n::Int = 5000, model::Union{Nothing, KaremakerModel} = nothing)
    rr = []
    symp_drive = []
    lf = []
    hf = []
    for T in Ts
        if model  === nothing
            model = KaremakerModel(T_resp = T)
        else
            model.T_resp = T
        end
        rr_T, symp_drive_T, lf_T, hf_T = lfhf2D(Dref, runoffs, model = model)
        push!(rr, rr_T)
            push!(lf, lf_T)
        push!(hf, hf_T)
        push!(symp_drive, symp_drive_T)
    end
    return rr, symp_drive, lf, hf
end