-
Notifications
You must be signed in to change notification settings - Fork 53
/
kernels.jl
100 lines (85 loc) · 4.03 KB
/
kernels.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
using GaussianProcesses, BenchmarkTools
using Random, Printf
Random.seed!(1)
const d = 10 # Input observation dimension
const nt = 100 # Number of training points
const λ = rand(d)
const X = 2π * rand(d, nt) # Training data
const cK = Array{Float64}(undef, nt, nt)
const kernelswithtags = Dict((k => ["simple", "isotropic"] for k in
(SEIso(1.0, 1.0), Mat12Iso(1.0, 1.0), Mat32Iso(1.0, 1.0),
Mat52Iso(1.0, 1.0), RQIso(1.0, 1.0, 1.0),
Periodic(1.0, 1.0, 1.0)))...,
(k => ["simple", "ARD"] for k in
(SEArd(λ, 1.0), Mat12Ard(λ, 1.0), Mat32Ard(λ, 1.0),
Mat52Ard(λ, 1.0), RQArd(λ, 1.0, 1.0)))...,
(k => ["simple", "non-stationary"] for k in
(LinIso(0.0), LinArd(λ)))...,
SEIso(1.0, 1.0) + RQIso(1.0, 1.0, 1.0) =>
["composite", "isotropic + isotropic"],
SEIso(1.0, 1.0) + LinIso(1.0) =>
["composite", "isotropic + non-stationary"],
Mat32Iso(1.0, 1.0) * Periodic(1.0, 1.0, 1.0) =>
["composite", "isotropic * isotropic"],
RQArd(λ, 1.0, 1.0) * LinIso(1.0) =>
["composite", "ARD * non-stationary"])
# SEIso(1.0, 1.0) + Noise(1.0) => ["composite", "isotropic + Noise"])
# Define benchmarks
const suite = BenchmarkGroup(["kernels"])
for (kernel, tags) in kernelswithtags
data = GaussianProcesses.KernelData(kernel, X, X)
np = GaussianProcesses.num_params(kernel)
stack = Array{Float64}(undef, nt, nt, np)
group = suite[string(typeof(kernel))] = BenchmarkGroup(tags)
group["data"] = @benchmarkable GaussianProcesses.KernelData($kernel, $X, $X)
group["cov"] = @benchmarkable GaussianProcesses.cov!($cK, $kernel, $X, $data)
group["grad"] = @benchmarkable GaussianProcesses.grad_stack!($stack, $kernel, $X, $data)
end
println("Benchmark suite:")
# TODO: Workaround until https://github.com/JuliaCI/BenchmarkTools.jl/issues/96 is fixed
show(stdout, MIME"text/plain"(), suite; verbose = true, limit = max(3, length(suite)))
println()
# If a cache of tuned parameters already exists, use it, otherwise, tune and cache
# the benchmark parameters. Reusing cached parameters is faster and more reliable
# than re-tuning `suite` every time the file is included.
const paramspath = joinpath(dirname(@__FILE__), "params.json")
if isfile(paramspath)
loadparams!(suite, BenchmarkTools.load(paramspath)[1], :evals, :samples)
else
println("\nTuning benchmarks...")
tune!(suite, verbose=true)
BenchmarkTools.save(paramspath, params(suite))
println()
end
println("\nRunning benchmarks...")
# Use @tagged to select more specific tests:
#const test_group = suite[@tagged "non-stationary"]
#const results = run(test_group, verbose = true, seconds = 1)
const results = run(suite, verbose = true, seconds = 1)
# Save results
println("\nSaving results...")
BenchmarkTools.save(joinpath(dirname(@__FILE__), "results.json"), results)
# Print table of results
function printpadright(text, total)
print(text)
for _ in (length(text)+1):total
print(" ")
end
end
function printresults(results)
# Compute max number of chars in the first column with kernel names
maxchar = max(6, maximum(length(k) for k in keys(results)))
# Print header
printpadright("kernel", maxchar)
println(" | data [μs] | cov [μs] | grad [μs]")
foreach(x -> print("-"), 1:maxchar)
println("-+-----------+-----------+----------")
# Print results for every kernel
foreach(results) do (key, result)
printpadright(key, maxchar)
@printf(" | %9.4f | %9.4f | %9.4f\n", time(result["data"]) / 1000,
time(result["cov"]) / 1000, time(result["grad"]) / 1000)
end
end
println("\nResults:")
printresults(results)