library(lpSolve)

maxim <- c(125,95)
LHS <- matrix(c(8,7,
                5,12)
              ,ncol=2, byrow=TRUE)
arah_kendala <- c("<=","<=")
RHS <- c(56,60)


# Penyelesaian
model <- lp ("max",
            maxim,
            LHS,
            arah_kendala,
            RHS,
            compute.sens = TRUE)

print(model)

model$solution

from scipy.optimize import linprog
obj = [-125,-95]

lhs_ineq = [[8,7],
           [5,12]]

rhs_ineq = [56,
           60]
bnd = [(0, float("inf")), 
       (0, float("inf"))]
model_2 = linprog(c=obj,
                 A_ub = lhs_ineq,
                 b_ub = rhs_ineq,
                 method = "revised simplex")

print(model_2)

print("nilai x dan y adalah = {}".format(model_2.x))

print("Untung Maksimumnya adalah = {}".format(model_2.fun*-1))

library(nloptr)
eval_f <- function( x ) {
return( list( "objective" = x[1]*x[4]*(x[1] + x[2] + x[3]) + x[3],
"gradient" = c( x[1] * x[4] + x[4] * (x[1] + x[2] + x[3]),
x[1] * x[4],
x[1] * x[4] + 1.0,
x[1] * (x[1] + x[2] + x[3]) ) ) )}

# Inequality constraints.
eval_g_ineq <- function( x ) {
constr <- c( 25 - x[1] * x[2] * x[3] * x[4] )
grad <- c( -x[2]*x[3]*x[4],
-x[1]*x[3]*x[4],
-x[1]*x[2]*x[4],
-x[1]*x[2]*x[3] )
return( list( "constraints"=constr, "jacobian"=grad ) )}

# Equality constraints.
eval_g_eq <- function( x ) {
constr <- c( x[1]^2 + x[2]^2 + x[3]^2 + x[4]^2 - 40 )
grad <- c( 2.0*x[1],
2.0*x[2],
2.0*x[3],
2.0*x[4] )
return( list( "constraints"=constr, "jacobian"=grad ) )}

# Initial values.
x0 <- c( 1, 5, 5, 1 )
# Lower and upper bounds of control.
lb <- c( 1, 1, 1, 1 )
ub <- c( 5, 5, 5, 5 )
# Optimal solution.
solution.opt <- c(1.00000000, 4.74299963, 3.82114998, 1.37940829)

# Set optimization options.
local_opts <- list( "algorithm" = "NLOPT_LD_MMA",
"xtol_rel" = 1.0e-7 )
opts <- list( "algorithm" = "NLOPT_LD_AUGLAG",
"xtol_rel" = 1.0e-7,
"maxeval" = 1000,
"local_opts" = local_opts,
"print_level" = 0 )

# Do optimization.
res <- nloptr( x0 = x0,
eval_f = eval_f,
lb = lb,
ub = ub,
eval_g_ineq = eval_g_ineq,
eval_g_eq = eval_g_eq,
opts = opts )
print(res)

res <- nloptr(
x0 = x0,
eval_f = eval_f,
lb = lb,
ub = ub,
eval_g_ineq = eval_g_ineq,
opts = opts )
print(res)

import numpy as np
import matplotlib.pyplot as plt

def f(x):
    return x1*x4*(x1 + x2 + x3) + x3

def obj_fun(x):
    return (x[0] - 0.5) ** 2 + 0.7 * x[0] * x[1] \
        + 1.2 * (x[1] + 0.7) ** 2

def gradient_fun(x):
    return np.array([2 * (x[0] - 0.5) + 0.7 * x[1], \
        0.7 * x[0] + 2 * 1.2 * (x[1] + 0.7)])

import numpy as np
from scipy.optimize import minimize, line_search, BFGS, SR1
from numdifftools import Gradient

class DescentAlgorithm:
    
    def __init__(self, fun, gradient=None, hess=None, nd={},
                 wolfe_c1=1e-4, wolfe_c2=0.1, x_tol=1e-6,
                 f_tol=1e-6, max_iter=50, save_history=False):
        
        self.fun = fun
        
        if gradient is None:
            self.gradient = Gradient(fun, **nd)
        else:
            self.gradient = gradient
        
        self.hess = hess
        self.wolfe_coefs = wolfe_c1, wolfe_c2
        self.x_tol = x_tol
        self.f_tol = f_tol
        self.max_iter = max_iter
        self.save_history = save_history
        self.history = []
    
    def optimize(self, x0, *args, **kwargs):
        
        x0 = np.atleast_1d(x0).astype(float)
        self.history = []
        xk = x0.copy()
        fk = self.fun(x0, *args, **kwargs)
        gradk = self.gradient(x0, *args, **kwargs)
        
        fc, gc = 1, 1
        
        pk = self.prepare_initial_step(xk, fk, gradk, *args,
                                       **kwargs)
        
        advance_x, advance_f, advance_max = True, True, True
        k = 0
        
        if self.save_history:
            self.history.append({"x":xk, "f":fk, "grad":gradk})
        
        while (advance_x or advance_f) and (k <= self.max_iter):
            
            alpha, fc_, gc_, fnew, fk, gradnew \
                = line_search(self.fun, self.gradient,
                              xk, sk, gradk, fk, args=args,
                              c1=self.wolfe_coefs[0],
                              c2=self.wolfe_coefs[1],
                              maxiter=15)
            
            if alpha is None:
                alpha = 1
                fnew = self.fun(xk + alpha * sk, *args, **kwargs)
                gradnew = self.gradient(xk + alpha * sk, *args,
                                        **kwargs)
            
            xnew = xk + alpha * pk
            fc = fc + fc_
            gc = gc + gc_
            
            if gradnew is None:
                gradnew = self.gradient(xnew, *args, **kwargs)
            
            advance_f = abs(fnew - fk) > self.f_tol
            advance_x = np.linalg.norm(xnew - xk) > self.x_tol
            
            xk, fk, gradk, pk = \
                self.prepare_next_step(xk, fk, gradk, pk,
                                       xnew, fnew, gradnew,
                                       *args, **kwargs)
            k = k + 1
            
            if self.save_history:
                self.history.append({"x":xk, "f":fk, "grad":gradk})
            
            if np.linalg.norm(pk) < np.sqrt(np.finfo(float).eps):
                self.message = 'Negligible step'
                self.success = True
                break
        
        if not (advance_x or advance_f):
            self.success = True
            self.message = 'Tolerance reached'
            
        elif k > self.max_iter:
            self.success = False
            self.message = 'Max iterations reached'
        
        self.x = xk
        self.f = fk
        self.grad = gradk
        self.fc = fc
        self.gc = gc
        self.result = {"x":xk, "f":fk, "grad":gradk, "iter":k,
                       "message":self.message,
                       "success":self.success}
    
    def prepare_next_step(self, xk, fk, gradk, pk, xnew, fnew, 
                          gradnew, *args, **kwargs):
        pass
    
    def prepare_initial_step(self, xk, fk, gradk, *args, **kwargs):
        pass

class SteepestDescent(DescentAlgorithm):

    def prepare_next_step(self, xk, fk, gradk, pk, xnew, fnew,
                          gradnew, *args, **kwargs):
        return xnew, fnew, gradnew, -gradnew

    def prepare_initial_step(self, xk, fk, gradk, *args, **kwargs):
        return -gradk

steepest = SteepestDescent(obj_fun, gradient=gradient_fun,
                           save_history=True)

library(tidyquant) 
library(plotly) 
library(timetk)

tick = c('BBCA.JK','TLKM.JK','BBRI.JK','ASII.JK','UNVR.JK')

pricedata = tq_get(tick,
                   from = '2020-01-01',
                   to   = '2022-12-09',              
                   get  = 'stock.prices')
pricedata

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import pandas
from pandas_datareader import data as pdr
import yfinance as yfin
yfin.pdr_override()
import datetime

tick1 = ['BBCA.JK','TLKM.JK','BBRI.JK','ASII.JK','UNVR.JK']
price_data=pdr.get_data_yahoo(tick1,
                              start='2020-01-01',
                              end= datetime.date.today())['Adj Close']

price_data

from numpy import log 
log_ret = np.log(price_data/price_data.shift(1))
log_ret

cov_mat = log_ret.cov() * 252
print(cov_mat)

# Simulating 5000 portfolios
num_port = 5000
# Creating an empty array to store portfolio weights
all_wts = np.zeros((num_port, len(price_data.columns)))
# Creating an empty array to store portfolio returns
port_returns = np.zeros((num_port))
# Creating an empty array to store portfolio risks
port_risk = np.zeros((num_port))
# Creating an empty array to store portfolio sharpe ratio
sharpe_ratio = np.zeros((num_port))

for i in range(num_port):
  wts = np.random.uniform(size = len(price_data.columns))
  wts = wts/np.sum(wts)
  
  # saving weights in the array
  
  all_wts[i,:] = wts
  
  # Portfolio Returns
  
  port_ret = np.sum(log_ret.mean() * wts)
  port_ret = (port_ret + 1) ** 252 - 1
  
  # Saving Portfolio returns
  
  port_returns[i] = port_ret
  
  
  # Portfolio Risk
  
  port_sd = np.sqrt(np.dot(wts.T, np.dot(cov_mat, wts)))
  
  port_risk[i] = port_sd
  
  # Portfolio Sharpe Ratio
  # Assuming 0% Risk Free Rate
  
  sr = port_ret / port_sd
  sharpe_ratio[i] = sr

names = price_data.columns
min_var = all_wts[port_risk.argmin()]
print(min_var)

max_sr = all_wts[sharpe_ratio.argmax()]
print(max_sr)

print(port_risk.min())

print(sharpe_ratio.max())

import plotly.express as px

# Variansi minimum
min_var = pd.Series(min_var, index=names)
min_var = min_var.sort_values()

fig = px.bar(min_var, color=min_var.index, title="Portofolio dengan Variansi minimum")
fig.show()

import plotly.express as px

# Maksimum Sharpe Ratio
max_sr = pd.Series(max_sr, index=names)
max_sr = max_sr.sort_values()

fig = px.bar(max_sr, color=max_sr.index, title="Portofolio dengan Tangensi")
fig.show()

fig = plt.figure()
ax1 = fig.add_axes([0.1,0.1,0.8,0.8])
ax1.set_xlabel('Risk')
ax1.set_ylabel("Returns")
ax1.set_title("Portfolio optimization and Efficient Frontier")
plt.scatter(port_risk, port_returns)
plt.show()