Final Exam Mat. Akt

Import Data

# Data
year <- c(1980, 1981, 1982, 1983, 1984, 1985, 1986, 1987, 1988, 1989,
          1990, 1991, 1992, 1993, 1994, 1995, 1996, 1997, 1998, 1999,
          2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009,
          2010, 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2019,
          2020, 2021, 2022, 2023)

productivity <- c(4169.442, 3951.818, 3502.381, 3750.342, 3528.127, 3344.306, 
                  3815.466, 4041.731, 3812.039, 3869.066, 3858.184, 3824.65, 
                  3853.345, 3824.334, 3947.745, 3921.171, 3961.818, 3661.605, 
                  4460.42, 4716.56, 4398.982, 4334.105, 4970.317, 5195.834, 
                  4964.834, 5101.404, 5212.358, 5220.077, 5389.627, 5391.495, 
                  5675.321, 5176.68, 5595.152, 5443.138, 5232.585, 5198.858, 
                  5325.21, 5796.57, 5819.193, 5863.358, 2305.98, 3483.08, 
                  6012.53, 5967.92)

temperature <- c(26, 26.6, 26.4, 26.7, 26.1, 26.3, 26.3, 26.6, 26.5, 26.4, 
                 26.6, 26.4, 26.5, 26.5, 26.4, 26.6, 26.5, 26.6, 27.1, 26.4, 
                 26.5, 26.7, 26.9, 26.8, 26.8, 26.9, 26.8, 26.8, 26.6, 27.0, 
                 27.1, 26.8, 26.9, 27.0, 27.0, 27.1, 27.4, 27.0, 27.1, 27.2, 
                 27.0, 27.3, 27.0, 27.1)

rainfall <- c(80.4, 78.3, 77.8, 78.6, 75.8, 78.4, 81.0, 79.7, 81.1, 79.3, 
              79.3, 79.8, 78.2, 76.6, 81.5, 78.2, 75.4, 76.4, 74.2, 80.4, 
              81.2, 82.5, 80.2, 75.9, 79.2, 81.0, 80.9, 81.2, 30.5, 80.3, 
              50.8, 60.2, 83.2, 50.5, 70.6, 70.9, 90.4, 80.6, 81.3, 90.9, 
              72.2, 70.6, 30.5, 80.3)

# Membuat dataframe
datafarm <- data.frame(year, productivity, temperature, rainfall)
print(datafarm)

##    year productivity temperature rainfall
## 1  1980     4169.442        26.0     80.4
## 2  1981     3951.818        26.6     78.3
## 3  1982     3502.381        26.4     77.8
## 4  1983     3750.342        26.7     78.6
## 5  1984     3528.127        26.1     75.8
## 6  1985     3344.306        26.3     78.4
## 7  1986     3815.466        26.3     81.0
## 8  1987     4041.731        26.6     79.7
## 9  1988     3812.039        26.5     81.1
## 10 1989     3869.066        26.4     79.3
## 11 1990     3858.184        26.6     79.3
## 12 1991     3824.650        26.4     79.8
## 13 1992     3853.345        26.5     78.2
## 14 1993     3824.334        26.5     76.6
## 15 1994     3947.745        26.4     81.5
## 16 1995     3921.171        26.6     78.2
## 17 1996     3961.818        26.5     75.4
## 18 1997     3661.605        26.6     76.4
## 19 1998     4460.420        27.1     74.2
## 20 1999     4716.560        26.4     80.4
## 21 2000     4398.982        26.5     81.2
## 22 2001     4334.105        26.7     82.5
## 23 2002     4970.317        26.9     80.2
## 24 2003     5195.834        26.8     75.9
## 25 2004     4964.834        26.8     79.2
## 26 2005     5101.404        26.9     81.0
## 27 2006     5212.358        26.8     80.9
## 28 2007     5220.077        26.8     81.2
## 29 2008     5389.627        26.6     30.5
## 30 2009     5391.495        27.0     80.3
## 31 2010     5675.321        27.1     50.8
## 32 2011     5176.680        26.8     60.2
## 33 2012     5595.152        26.9     83.2
## 34 2013     5443.138        27.0     50.5
## 35 2014     5232.585        27.0     70.6
## 36 2015     5198.858        27.1     70.9
## 37 2016     5325.210        27.4     90.4
## 38 2017     5796.570        27.0     80.6
## 39 2018     5819.193        27.1     81.3
## 40 2019     5863.358        27.2     90.9
## 41 2020     2305.980        27.0     72.2
## 42 2021     3483.080        27.3     70.6
## 43 2022     6012.530        27.0     30.5
## 44 2023     5967.920        27.1     80.3

Nomor 1

library(SciViews)
lnproductivity = ln(datafarm$productivity)
lntemperature = ln(datafarm$temperature)
lnrainfall = ln(datafarm$rainfall)

datafarmnew = data.frame(year, lnproductivity, lntemperature, lnrainfall)
datafarmnew

##    year lnproductivity lntemperature lnrainfall
## 1  1980       8.335537      3.258097   4.387014
## 2  1981       8.281931      3.280911   4.360548
## 3  1982       8.161198      3.273364   4.354141
## 4  1983       8.229602      3.284664   4.364372
## 5  1984       8.168522      3.261935   4.328098
## 6  1985       8.115014      3.269569   4.361824
## 7  1986       8.246818      3.269569   4.394449
## 8  1987       8.304428      3.280911   4.378270
## 9  1988       8.245919      3.277145   4.395683
## 10 1989       8.260768      3.273364   4.373238
## 11 1990       8.257952      3.280911   4.373238
## 12 1991       8.249222      3.273364   4.379524
## 13 1992       8.256697      3.277145   4.359270
## 14 1993       8.249140      3.277145   4.338597
## 15 1994       8.280900      3.273364   4.400603
## 16 1995       8.274146      3.280911   4.359270
## 17 1996       8.284458      3.277145   4.322807
## 18 1997       8.205657      3.280911   4.335983
## 19 1998       8.402998      3.299534   4.306764
## 20 1999       8.458835      3.273364   4.387014
## 21 2000       8.389128      3.277145   4.396915
## 22 2001       8.374270      3.284664   4.412798
## 23 2002       8.511239      3.292126   4.384524
## 24 2003       8.555612      3.288402   4.329417
## 25 2004       8.510135      3.288402   4.371976
## 26 2005       8.537271      3.292126   4.394449
## 27 2006       8.558788      3.288402   4.393214
## 28 2007       8.560267      3.288402   4.396915
## 29 2008       8.592231      3.280911   3.417727
## 30 2009       8.592578      3.295837   4.385770
## 31 2010       8.643882      3.299534   3.927896
## 32 2011       8.551919      3.288402   4.097672
## 33 2012       8.629656      3.292126   4.421247
## 34 2013       8.602111      3.295837   3.921973
## 35 2014       8.562661      3.295837   4.257030
## 36 2015       8.556194      3.299534   4.261270
## 37 2016       8.580207      3.310543   4.504244
## 38 2017       8.665022      3.295837   4.389499
## 39 2018       8.668917      3.299534   4.398146
## 40 2019       8.676478      3.303217   4.509760
## 41 2020       7.743261      3.295837   4.279440
## 42 2021       8.155672      3.306887   4.257030
## 43 2022       8.701601      3.295837   3.417727
## 44 2023       8.694154      3.299534   4.385770

Nomor 2

library(SciViews)
modelfarm <- lm(ln(productivity) ~ ln(rainfall) + ln(temperature), data = datafarm)

summary(modelfarm)

## 
## Call:
## lm(formula = ln(productivity) ~ ln(rainfall) + ln(temperature), 
##     data = datafarm)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.75263 -0.06896  0.01269  0.09988  0.19079 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -19.1642     7.4161  -2.584 0.013422 *  
## ln(rainfall)     -0.1766     0.1194  -1.479 0.146726    
## ln(temperature)   8.6217     2.2311   3.864 0.000389 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1744 on 41 degrees of freedom
## Multiple R-squared:  0.3161, Adjusted R-squared:  0.2827 
## F-statistic: 9.473 on 2 and 41 DF,  p-value: 0.0004149

persamaan nya adalah \[ ln(Y) = -19.164 + 8.622 ln(Temperature) - 0.1766 ln(Rainfall) \]

Dengan nilai R-squared sebesar 28.27% namun secara parameter sudah signifikan

# Menghitung rata-rata suhu
rata_rata_suhu <- mean(datafarmnew$lntemperature)

# Menghitung rata-rata curah hujan
rata_rata_curah_hujan <- mean(datafarmnew$lnrainfall)

# Menampilkan hasil
print(paste("Rata-rata suhu:", rata_rata_suhu))

## [1] "Rata-rata suhu: 3.28586892485743"

print(paste("Rata-rata curah hujan:", rata_rata_curah_hujan))

## [1] "Rata-rata curah hujan: 4.29938900569327"

E[G]

print(paste("𝝁𝑮 =", "𝒂+ 𝜷𝟏𝝁𝟏 + 𝜷𝟐𝝁𝟐"))

## [1] "𝝁𝑮 = 𝒂+ 𝜷𝟏𝝁𝟏 + 𝜷𝟐𝝁𝟐"

print(paste("𝝁𝑮 =", " -19.164 +  8.622 (3.2858) - 0.1766(4.299)"))

## [1] "𝝁𝑮 =  -19.164 +  8.622 (3.2858) - 0.1766(4.299)"

a = -19.164 + 8.622*3.2858 - 0.176*4.299
a

## [1] 8.409544

exp(a)

## [1] 4489.711

 mean(predict(modelfarm, type = "response"))

## [1] 8.406432

Ekspektasi Produksi menggunakan Model Linear adalah 4485.287

Nomor 3

# Menghitung variansi dari kolom Suhu
variansi_suhu <- var(datafarm$temperature)

# Menghitung variansi dari kolom Curah_Hujan
variansi_curah_hujan <- var(datafarm$rainfall)

variansi_gabungan = cov(datafarm$temperature, datafarm$rainfall)


# Menampilkan hasil
print(paste("Variansi Suhu:", variansi_suhu))

## [1] "Variansi Suhu: 0.103229386892178"

print(paste("Variansi Curah Hujan:", variansi_curah_hujan))

## [1] "Variansi Curah Hujan: 155.51971987315"

print(paste('Variansi Gabungan:', variansi_gabungan))

## [1] "Variansi Gabungan: -0.523134249471463"

Variance Model

variance_value = (((8.622)^2 *0.1032) + ((-0.1766)^2 * 155.52) + (2*(8.622*(-0.1766))*(-0.523)) + 0.01)
                                                            
    variance_value                                                      

## [1] 14.12475

Nomor 4

Premi Loading Faktor

Premi1 <- (1+0.01)*exp(a)
Premi2 <- (1+0.02)*exp(a)
Premi3 <- (1+0.03)*exp(a)
Premi4 <- (1+0.04)*exp(a)
Premi5 <- (1+0.05)*exp(a)
Premi6 <- (1+0.06)*exp(a)
Premi7 <- (1+0.07)*exp(a)
Premi8 <- (1+0.08)*exp(a)
Premi9 <- (1+0.09)*exp(a)
Premi10 <- (1+0.1)*exp(a)

Premi <- c(Premi1, Premi2, Premi3, Premi4, Premi5, Premi6, Premi7, Premi8, Premi9, Premi10)
LF <- c(0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1)
Datax <- data.frame(Loading_Faktor = LF, 
                    Premi = Premi)
Datax

##    Loading_Faktor    Premi
## 1            0.01 4534.608
## 2            0.02 4579.505
## 3            0.03 4624.402
## 4            0.04 4669.299
## 5            0.05 4714.196
## 6            0.06 4759.094
## 7            0.07 4803.991
## 8            0.08 4848.888
## 9            0.09 4893.785
## 10           0.10 4938.682

Premi Standar Deviasi

  UP <- 1
Premi11 <- ((0.01)* sqrt((variance_value)) +exp(a))*UP
Premi22 <- ((0.02)* sqrt((variance_value)) +exp(a))*UP
Premi33 <- ((0.03)* sqrt((variance_value)) +exp(a))*UP
Premi44 <- ((0.04)* sqrt((variance_value)) +exp(a))*UP
Premi55 <- ((0.05)* sqrt((variance_value)) +exp(a))*UP
Premi66 <- ((0.06)* sqrt((variance_value)) +exp(a))*UP
Premi77 <- ((0.07)* sqrt((variance_value)) +exp(a))*UP
Premi88 <- ((0.08)* sqrt((variance_value)) +exp(a))*UP
Premi99 <- ((0.09)* sqrt((variance_value))  +exp(a))*UP
Premi1010  <-((0.1)* sqrt((variance_value)) + exp(a))*UP



Datax$Premi_StDev = c(Premi11, Premi22, Premi33, Premi44, Premi55, Premi66, Premi77, Premi88, Premi99, Premi1010)
  
Datax

##    Loading_Faktor    Premi Premi_StDev
## 1            0.01 4534.608    4489.749
## 2            0.02 4579.505    4489.786
## 3            0.03 4624.402    4489.824
## 4            0.04 4669.299    4489.861
## 5            0.05 4714.196    4489.899
## 6            0.06 4759.094    4489.936
## 7            0.07 4803.991    4489.974
## 8            0.08 4848.888    4490.012
## 9            0.09 4893.785    4490.049
## 10           0.10 4938.682    4490.087

library(plotly)

## Loading required package: ggplot2

## 
## Attaching package: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

 plot_ly(x = Datax$Loading_Faktor, y = Datax$Premi, type = 'scatter', mode = 'lines', name = 'Premi LF') %>%
  add_trace( x = Datax$Loading_Faktor , y = Datax$Premi_StDev , mode = 'lines', name = 'Premi StaDev') %>%
  layout(title = 'Actual vs. ARIMA Forecast',
         xaxis = list(title = 'Loading Faktor'),
         yaxis = list(title = 'Premi'))

3D Plot Model Cobb Douglas

curah_hujan = exp(datafarmnew$lnrainfall)
suhu = exp(datafarmnew$lntemperature)
Produktivitas = exp((-19.164 + 8.622*(datafarmnew$lntemperature) - 0.176*(lntemperature)))
  
  
  
fig <- plot_ly(x = ~curah_hujan, y = ~suhu, z = ~Produktivitas , type = 'scatter3d', mode = 'lines',name = 'Coub Douglass Model' , line = list(width = 4, color = ~c, colorscale = list(c(0,'#FA52EB'), c(1,'#ADB240'))))  %>% 
  layout(title = '3d Plot on Coub-Douglas Model',
         xaxis = list(title = 'ln Curah Hujan'),
         yaxis = list(title = 'ln Suhu'), 
         zaxis = list(title = 'ln Produktivitas' ))

fig

## Warning: 'layout' objects don't have these attributes: 'zaxis'
## Valid attributes include:
## '_deprecated', 'activeshape', 'annotations', 'autosize', 'autotypenumbers', 'calendar', 'clickmode', 'coloraxis', 'colorscale', 'colorway', 'computed', 'datarevision', 'dragmode', 'editrevision', 'editType', 'font', 'geo', 'grid', 'height', 'hidesources', 'hoverdistance', 'hoverlabel', 'hovermode', 'images', 'legend', 'mapbox', 'margin', 'meta', 'metasrc', 'modebar', 'newshape', 'paper_bgcolor', 'plot_bgcolor', 'polar', 'scene', 'selectdirection', 'selectionrevision', 'separators', 'shapes', 'showlegend', 'sliders', 'smith', 'spikedistance', 'template', 'ternary', 'title', 'transition', 'uirevision', 'uniformtext', 'updatemenus', 'width', 'xaxis', 'yaxis', 'barmode', 'bargap', 'mapType'