#Problem 1 #definition of a straight line function linrel=function(x,m,b) m*x+b #definition x-Vctor x_vec=seq(-1,3,0.1) #definition y-vector y_vec=linrel(x_vec,1,0) #plotting x-vector vs. y-vector plot(x_vec, y_vec, type="l", col="black", xlab='x',ylab='f(x)=y') #generating mock data point for 0 <= x <= 2 for 10 equaly spaced points equ_points=10 mock_vec = rep(NA,10) x_vec_2 =seq(0,2, 2/(equ_points-1)) mock_vec = rnorm({1:10},mean = x_vec_2, sd = 0.4) points(x_vec_2,mock_vec, col="red", type="p") ##Problem2 conrel =function(b) b quarel =function(x,a,m,b) a*x**2 + m*x+b reso = 50 #the resolution a_vec = seq(-4,4,length.out = reso) m_vec = seq(-4,8,length.out = reso) b_vec = seq(-2,2,length.out = reso) # Define function chi2_quad chi2_qua = function(a,m,b) sum((mock_vec - quarel(mock_vec,a,m,b))**2/(0.4)**2) chi2_lin = function(m,b) sum((mock_vec - linrel(mock_vec,m,b))**2/(0.4)**2) chi2_con = function(b) sum((mock_vec - conrel(b))**2/(0.4)**2) # Calculate 3-dim grid for values a=-4..4, m=-4..8, b-2..2 qua_arr = array(NA, dim=c(length(b_vec),length(m_vec),length(a_vec))) lin_arr = array(NA, dim=c(length(b_vec),length(m_vec))) con_arr = array(NA, dim=c(length(b_vec))) for (i in (1:length(b_vec))){ con_arr[i] = chi2_con(b_vec[i]) for (j in (1:length(m_vec))){ lin_arr[i,j] = chi2_lin(m_vec[j],b_vec[i]) for (k in (1:length(a_vec))){ qua_arr[i,j,k] = chi2_qua(a_vec[k],m_vec[j],b_vec[i]) } } } #find min minimum_qua=min(qua_arr) best_m_qua = m_vec[which(qua_arr == minimum_qua,arr.ind=TRUE)[2]] best_b_qua = b_vec[which(qua_arr == minimum_qua,arr.ind=TRUE)[1]] best_a_qua = a_vec[which(qua_arr == minimum_qua,arr.ind=TRUE)[3]] points(x_vec, quarel(x_vec,best_a_qua,best_m_qua,best_b_qua), type="l", col="red", xlab='x',ylab='f(x)=y') minimum_lin=min(lin_arr) best_m_lin = m_vec[which(lin_arr == minimum_lin,arr.ind=TRUE)[2]] best_b_lin = b_vec[which(lin_arr == minimum_lin,arr.ind=TRUE)[1]] points(x_vec, linrel(x_vec,best_m_lin,best_b_lin), type="l", col="blue", xlab='x',ylab='f(x)=y') minimum_con=min(lin_arr) best_b_con = b_vec[which(con_arr == minimum_con,arr.ind=TRUE)[1]] points(x_vec, rep(conrel(best_b_lin),41), type="l", col="yellow", xlab='x',ylab='f(x)=y') #Probelm3 norm = function(a) 1/sum(a) * a posterior_qua = exp(-0.5 *(norm(qua_arr) - minimum_qua)) posterior_lin = exp(-0.5 *(norm(lin_arr) - minimum_lin)) posterior_con = exp(-0.5 *(norm(con_arr) - minimum_con))