chapter9.jl
# Chapter 9
#import Pkg
#Pkg.add("Plots")
#Pkg.add("Optim")
#Pkg.add("StatsKit")
#Pkg.add("GLM")
#Pkg.add("DataFrames")
#Pkg.add("StatsPlots")
#Pkg.add("StatFiles")
#Pkg.add("NLSolversBase")
#Pkg.add("Distributions")
using Random
using Statistics
using Plots
using Optim
using NLSolversBase
using StatsKit
using CSV
using GLM
using DataFrames
using LinearAlgebra
using StatsPlots
using StatFiles
using Distributions
# 9.1
# 9.2
# 9.2.1
# 9.2.2
Random.seed!(123456789)
rng = MersenneTwister(123456789)
M = 1000
data1 = Array{Union{Missing, Float64}}(missing, M, 12)
for i = 1:M
N = rand(2:10, 1)[1]
v = rand!(rng, zeros(N))
b = (N - 1).*v./N
p = maximum(b)
x = Vector{Union{Missing, Float64}}(missing, 10)
x[1:N] = b
data1[i,1] = N
data1[i,2] = p
data1[i,3:12] = x
#println(i)
end
data1 = DataFrame(Num = data1[1:end,1], Price = data1[1:end,2], Bid1 = data1[1:end,3], Bid2 = data1[1:end,4], Bid3 = data1[1:end,5], Bid4 = data1[1:end,6], Bid5 = data1[1:end,7], Bid6 = data1[1:end,8], Bid7 = data1[1:end,9], Bid8 = data1[1:end,10], Bid9 = data1[1:end,11], Bid10 = data1[1:end,12])
# 9.2.3
# 9.2.4
sealed_2bid = function(bids; ϵ = 0.5)
values = zeros(length(bids))
for i = 1:length(bids)
Ĥ = mean(bids .< bids[i])
ĥ = (mean(bids .< (bids[i] + ϵ)) - mean(bids .< (bids[i] - ϵ)))./(2*ϵ)
values[i] = bids[i] + Ĥ./ĥ
end
return values
end
bids = [data1[data1.Num .== 2,3];data1[data1.Num .== 2,4]]
values = sealed_2bid(bids)
density(bids)
density!(values)
hline!(fill(0.5, 206))
xlabel!("value")
ylabel!("density")
# 9.3
# 9.3.1
# 9.3.2
# 9.3.3
# 9.3.4
English_2bid = function(price; K=100, ϵ = 1e-8)
min1 = minimum(price)
max1 = maximum(price)
diff1 = (max1 - min1)/K
Fᵥ = Array{Union{Missing, Float64}}(missing, K, 2)
min_temp = min1 - ϵ
max_temp = min_temp + diff1
Fᵥ[1,1] = (min_temp + max_temp)/2
gp = mean((price .> min_temp).*(price .< max_temp))
Fᵥ[1,2] = gp/2
for k = 2:K
min_temp = max_temp - ϵ
max_temp = min_temp + diff1
Fᵥ[k,1] = (min_temp + max_temp)/2
gp = mean((price .> min_temp).*(price .< max_temp))
Fᵥ[k,2] = gp/(2*(1 - Fᵥ[k-1,2])) + Fᵥ[k-1,2]
end
return Fᵥ
end
Random.seed!(123456789)
rng = MersenneTwister(123456789)
M = 10000
data2 = Array{Union{Missing, Float64}}(missing, M, 2)
for i = 1:M
N = rand(2:10, 1)[1]
v = randn!(rng, zeros(N))
b = v
p = -sort(-b)[2]
data2[i,1] = N
data2[i,2] = p
end
data2 = DataFrame(Num = data2[1:end,1], Price = data2[1:end,2])
price = data2.Price[data2.Num .== 2]
Fᵥ = English_2bid(price)
plot(Fᵥ[1:end,1],Fᵥ[1:end,2])
plot!(Fᵥ[1:end,1],cdf(Normal(),Fᵥ[1:end,1]))
# 9.4
# 9.4.1
path = "ENTER PATH HERE"
x = DataFrame(CSV.read(string(path,"/auctions.csv")))
log_value = Vector{Union{Missing, Float64}}(missing, 8920)
for i = 1:8920
if !(x.log_value[i] == "NA")
log_value[i] = parse(Float64, x.log_value[i])
end
end
x.log_value = log_value
x.Sale_Size = x[1:end,11]
x.Road_Construction = x[1:end,13]
x.Method_of_Sale = x[1:end,10]
x.Species = string.(x.Species)
x.Region = string.(x.Region)
x.Forest = string.(x.Forest)
x.District = string.(x.District)
x1 = dropmissing(x)
lm1 = lm(@formula(log_amount ~ Salvage + Acres + log_value + Sale_Size + Road_Construction + Species + Region + Forest + District), x1)
x1.norm_bid = residuals(lm1)
histogram(x1.norm_bid)
xlabel!("bid")
ylabel!("Frequency")
# 9.4.2
y = x1[(x1.num_bidders .== 2).*(x1.Method_of_Sale .== "S"),1:end]
bids = y.norm_bid
values = sealed_2bid(bids)
quantile(bids, [0, 0.25, 0.5, 0.75, 1])
mean(bids)
quantile(values, [0, 0.25, 0.5, 0.75, 1])
mean(values)
# 9.4.3
y = x1[(x1.num_bidders .== 2).*(x1.Method_of_Sale .== "A"),1:end]
price = y.norm_bid[y.Rank .== 2]
Fᵥ = English_2bid(price)
# 9.4.4
scatter(Fᵥ[1:end,1],Fᵥ[1:end,2])
scatter!(ecdf(values[values .< 1.5]))
# 9.5
# 9.5.1
scatter(x1.num_bidders[x1.Method_of_Sale .== "A"],x1.norm_bid[(x1.Method_of_Sale .== "A").*(x1.Rank .== 2)])
# 9.5.2
# 9.5.3
# 9.5.4
English_Nbid = function(price, N; K=100, ϵ = 1e-8)
min1 = minimum(price)
max1 = maximum(price)
diff1 = (max1 - min1)/K
Fᵥ = Array{Union{Missing, Float64}}(missing, K, 2)
min_temp = min1 - ϵ
max_temp = min_temp + diff1
Fᵥ[1,1] = (min_temp + max_temp)/2
gp = mean((price .> min_temp).*(price .< max_temp))
Fᵥ[1,2] = (gp/N)^(1/(N-1))
for k = 2:K
min_temp = max_temp - ϵ
max_temp = min_temp + diff1
Fᵥ[k,1] = (min_temp + max_temp)/2
gp = mean((price .> min_temp).*(price .< max_temp))
Fᵥ[k,2] = gp/(N*(N-1)*((Fᵥ[k-1,2]^(N-2)))*(1 - Fᵥ[k-1,2])) + Fᵥ[k-1,2]
end
return Fᵥ
end
y₆ = x1[(x1.num_bidders .> 5).*(x1.Method_of_Sale .== "A"),1:end]
# English auctions with at least 6 bidders.
price_6 = y₆.norm_bid[y₆.Rank .== 2]
Fᵥ_6 = English_Nbid(price_6, 6)
Fᵥ_3 = English_Nbid(price_6, 3)
Fᵥ_2 = English_2bid(price_6)
# 9.5.5
scatter(Fᵥ[1:end,1],Fᵥ[1:end,2], label = "Values (2 bidders)")
scatter!(Fᵥ_6[1:end,1],Fᵥ_6[1:end,2], label = "Values (6 bidders)")
scatter!(Fᵥ_3[1:end,1],Fᵥ_3[1:end,2], label = "Values (as if 3 bidders)")
scatter!(Fᵥ_2[1:end,1],Fᵥ_2[1:end,2], label = "Values (as if 2 bidders)", legend = :topleft)
# 9.5.6
# 9.5.7
sealed_Nbid = function(bids, N; ϵ = 0.5)
values = zeros(length(bids))
for i = 1:length(bids)
Ĥ = mean(bids .< bids[i])
ĥ = (mean(bids .< (bids[i] + ϵ)) - mean(bids .< (bids[i] - ϵ)))./(2*ϵ)
Ĝ = Ĥ.^(N-1)
ĝ = (N-1).*ĥ.*(Ĥ.^(N-2))
values[i] = bids[i] + Ĝ./ĝ
end
return values
end
y₆ = x1[(x1.num_bidders .> 5).*(x1.Method_of_Sale .== "S"), 1:end]
bids_6 = y₆.norm_bid
values_6 = sealed_Nbid(bids_6, 6, ϵ=3)
scatter(Fᵥ[1:end,1],Fᵥ[1:end,2], label = "Values (2 bidder English)")
scatter!(ecdf(values[values .< 1.5]), label = "Values (2 bidder sealed)")
scatter!(ecdf(values_6[values_6 .< 1.5]), label = "Values (6 bidder sealed)", legend = :topleft)
# 9.6