Open shreyam202236 opened 1 year ago
knitr::opts_chunk$set(echo = TRUE)
knitr::opts_chunk$set(echo = TRUE) knitr::opts_chunk$set(fig.align = 'center') knitr::opts_chunk$set(fig = TRUE) knitr::opts_chunk$set(warning = FALSE) knitr::opts_chunk$set(tidy.opts = list(width.cutoff = 60), tidy = TRUE)
library(tidyverse)
gaussian <- function(vec) { b0 <- vec[1] b1 <- vec[2] sigma <- vec[3] loglike <- 0 for (i in 1:n row(Insurance)) { loglike <- loglike + dnorm(Insurance$Claims[i],mean=avg, sd=sigma, log=TRUE) avg <- b0 + b1*Insurance$Holders[i] } return(-1*loglike) } bic <- function(mle, parameters, data) { N <- length(data) return((length(parameters))log(N)-2(mle$value)) } gauss.est <- optim(c(0,1/5,1), gaussian) gauss.bic <- bic(gauss.est, gauss.est$par, Insurance) cat("The Gaussian BIC is:", gauss.bic)
loglap <- function(x, mu, b) { deno <- 2*b num <- exp(-abs(x-mu)/b) laplace <- num/deno return(log(laplace)) } laplacian <- function(vec) { b0 <- vec[1] b1 <- vec[2] sigma <- vec[3] loglike <- 0 for (i in 1:nrow(Insurance)) { avg <- b0 + b1*Insurance$Holders[i] loglike <- loglike + loglap(Insurance$Claims[i],mu=avg, b=sigma) } return(-1*loglike) } laplac.est <- optim(c(0, 0, 1), laplacian) laplace.bic <- bic(laplac.est, laplac.est$par, Insurance) cat("The Laplacian BIC is:", laplac.bic)
lognorm <- function(vec) { b0 <- vec[1] b1 <- vec[2] sigma <- vec[3] loglike <- 0 for (i in 1:nrow(Insurance)) { if (Insurance$Claims[i] > 0) { avg <- b0 + b1*Insurance$Holders[i] loglike <- loglike + dlnorm(Insurance$Claims[i], meanlog=avg, sdlog=sigma, log=TRUE) } } return(-1*loglike) } lognorm.est <- optim(c(0, 0, 1), lognormal) lognorm.bic <- bic(lognorm.est, lognorm.est$par, Insurance) cat("The Lognormal BIC is:", lognorm.bic)
gamma.mod <- function(vec) { b0 <- vec[1] b1 <- vec[2] sigma <- vec[3] loglike <- 0 for (i in 1:nrow(Insurance)) { if (Insurance$Claims[i] > 0) { avg <- b0 + b1*Insurance$Holders[i] loglike <- loglike + dgamma(Insurance$Claims[i], shape=avg, scale=sigma, log=TRUE) } } return(-1*loglike) } gamma.est <- optim(c(0,1/5,1), gamma.mod) gamma.bic <- bic(gamma.est, gamma.est$par, Insurance) cat("The Gamma BIC is:", gamma.bic)
Checked but there were no review by other ream member
title: "Problem 4" author: "Dipanjoy Saha" date: "2022-11-17" output: html_document
Problem 4