def sim_prediction(self,mu,Y,h,t_params,simulations):		

		sim_vector = np.zeros([simulations,h])

		for n in range(simulations):
			# Create arrays to iteratre over		
			Y_exp = copy.deepcopy(Y)
			mu_exp = copy.deepcopy(mu)

			# Loop over h time periods			
			for t in range(h):

				new_value = t_params[0] + np.random.randn(1)*t_params[len(t_params)-1]

				if self.ar != 0:
					for j in range(1,self.ar+1):
						new_value += t_params[j]*Y_exp[len(Y_exp)-j]

				if self.ma != 0:
					for k in range(1,self.ma+1):
						if (k-1) >= t:
							new_value += t_params[k+self.ar]*(Y_exp[len(Y_exp)-k]-mu_exp[len(mu_exp)-k])

				Y_exp = np.append(Y_exp,[new_value])
				mu_exp = np.append(mu_exp,[0]) # For indexing consistency

				sim_vector[n] = Y_exp[(len(Y_exp)-h):(len(Y_exp))]

		return np.transpose(sim_vector)
Simulates a h-step ahead mean prediction
dataset/ETHPy150Open RJT1990/pyflux/pyflux/arma/arma.py/ARIMA.sim_prediction