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