Abstrak —PT Japfa Comdeed Indonesia Tbk merupakan salah satu  perusahaan  yang  bergerak  di  bidang  agrofood.  PT  JCI Tbk  mengolah  limbah  air  hasil  produksinya  berdasarkan Peraturan  Gubernur  Jawa  Timur  No.  72  Tahun  2013  tentang baku  mutu  air  limbah  bagi  industri  dan/atau  kegiatan  usaha lainnya.  Meskipun  telah  dilakukan  pengolahan  limbah  yang sesuai namun masih ditemukan pengamatan pengolahan limbah yang  melebihi  dari  baku  mutu  yang  telah  ditetapkan.  Selama ini PT JCI Tbk hanya melakukan pencacatan dalam memonitor proses pengolahan limbah, sehingga perlu dilakukan pengendalian  kualitas  pengolahan  limbah  secara  statistik. Penelitian ini menggunakan dua variabel yaitu Total Dissolved Solid  (TDS)  dan  pH.  Adanya  autokorelasi  menyebabkan  false alarm  pada  grafik  kendali  sehingga  pada  penelitian  ini digunakan  grafik  kendali  multivariat  bersasis  time  series. Metode  yang  digunakan  digunakan  adalah  grafik  kendali kombinasi  MEWMA  dimana  grafik  kendali  M 2 Z 2 digunakan untuk  memonitor  variabilitas  proses  dan  peta  kenali  MZ digunakan  untuk  memonitor  mean  proses.  Pada  analisis  time series diperoleh model terbaik yang terpilih  adalah VAR([1,5,8])  sehingga  residual  dari  model tersebut digunakan untuk mengendalikan proses pengolahan limbah air. Pada fase I  variabilitas  proses  dan  mean  proses  telah  stabil,  dan  pada fase  II  variabilitas  proses  telah  stabil  namun  mean  proses belum stabil. Kapabilitas proses diperoleh nilai MPp 1,405 dan MPpk  0,875  sehingga  dapat  dikatakan  bahwa  proses pengolahan  limbah  memiliki  presisi  yang  baik  namun akurasinya rendah.
Kata Kunci— Autokorelasi, Grafik Kendali Kombinasi MEWMA, Limbah Air, pH, TDS

## I. PENDAHULUAN
T. Japfa Comfeed Indonesia (JCI) Tbk  merupakan salah satu perusahaan besar yang bergerak dibidang agro- food dengan  hasil  produksi  berupa  pakan  ternak  untuk  ayam  ras pedaging.  Ayam  ras  pedaging  merupakan  salah  satu komoditas unggas yang peminatnya tinggi, dengan tingginya permintaan akan ayam menyebabkan meningkatnya permintaan pakan ternak. Tingginya permintaan pakan ternak membuat  PT  JCI  Tbk  berusaha  agar  dapat  memenuhi permintaan  dengan  meningkatkan  aktivitas  produksi. Semakin  tinggi  produksi  yang  dihasilkan  juga  berdampak pada  banyaknya  limbah  yang  dihasilkan.  PT  JCI  Tbk mengolah limbah air hasil produksinya berdasarkan Peraturan Gubernur  Jawa  Timur  No.  72  Tahun  2013  tentang  baku mutu  air  limbah  bagi  industri  dan/atau  kegiatan  usaha lainnya.  Limbah  hasil  dari  proses  produksi  sebelum  dibuang ke  saluran  harus  diolah  hingga  limbah  aman  dan  tidak merusak  lingkungan.  Selama  ini  PT  JCI  Tbk  melakukan pengontrolan  terhadap  limbah  air  secara  deskriptif  untuk
menentukan  apakah  limbah  layak  untuk  dibuang  ke  saluran atau  diperlukan  proses  kembali  pengolahan  limbah  ketika nilai pH dan Total Dissolved Solid (TDS) melebihi dari baku mutu yang ditentukan. Sehingga pada penelitian ini saya akan menggunakan  metode  statistika  yaitu  grafik  kendali multivariat untuk pengendalian kualitas limbah air.
Penelitian ini akan menggunakan grafik kendali kombinasi MEWMA  dengan  pendekatan time  series karena  pemberian larutan  kimiawi  yang  dilakukan  memiliki  pengaruh  antara pengamatan  saat  ini  dengan  pengamatan  sebelumnya sehingga  adanya  indikasi  autokorelasi  pada  data  kualitas pengolahan  limbah  yang  akan  diteliti.  Reynolds  &  Cho  [1] dalam  penelitiannya  yang  berjudul “Multivariate  Control Chart  for  Monitoring  The  Mean  Vector  and  Covariance Matrix” menyebutkan  bahwa  grafik  kendali  kombinasi MEWMA  memiliki  hasil  yang  lebih  baik  jika  dibandingkan dengan  grafik  kendali  MEWMA  standart  dan  grafik  kendali T 2 Hotteling  berdasarkan Average  Time  to  Signal (ATS). Grafik  kendali  kombinasi  MEWMA  menggunakan  grafik kendali  M 2 Z 2 untuk  memonitor variabiltas proses dan grafik kendali  MZ  untuk  memonitor mean proses.  Data  yang digunakan  merupakan  data  pengolahan  limbah  air  mulai Januari  2017  –  Februari  2019.  Karakteristik  kualitas  yang digunakan adalah TDS dan pH.
Rumusan  masalah  pada  penelitian  ini  adalah  bagaimana karakteristik  data  dari  karakteristik  kualitas  limbah  air  PT JCI,  bagaimana  model  terbaik  yang  terbentuk  dari  analisis time  series dan  bagaimana  kapabilitas  proses  pengolahan limbah  air  PT  JCI.  Tujuan  penelitian  ini  menjawab  dari rumusan  masalah.  Manfaat yang diperoleh dari penelitian ini diharapkan  dapat  menjadi  bahas  evaluasi  PT  Japfa  Comfeed Indonesia  Tbk  Margomulyo  dalam  peningkatan  pengolahan limbah industri.

## II. TINJAUAN PUSTAKA
A. Model Vector Autoregressive (VAR)
Model Vector  Autoregressive (VAR)  merupakan  model hasil  pengembangan  dari Autoregressive (AR)  yang  dapat diidentifikasi  menggunakan  plot Partial  Cross-Correlation Matrix Function. Model VAR digunakan untuk memodelkan dua  variabel  atau  lebih  yang  memiliki  hubungan  satu  sama lain.
t 1 t-1 p t-p t Z = μ + Φ Z + …+ Φ Z + a

## B. Stasioneritas Data
Stasioneritas  dalam  data time  series ditunjukkan  apabila rata-rata  dan  varians berfluktuasi konstan setiap waktu. Jika plot time  series berfluktuasi  di  sekitar  garis  yang  sejajar
Pengendalian Kualitas Pengolahan Limbah Air PT. Japfa Comfeed Indonesia Tbk Menggunakan Grafik Kendali

## Multivariat Berbasis Time Series
Ilma Tamarina Arba, Wibawati, Ni Luh Putu Satyaning P. P. Departemen Statistika, Fakultas Matematika, Komputasi, dan Sains Data, Institut Teknologi Sepuluh Nopember (ITS)
e-mail : wibawati.its@gmail.com
P
dengan  sumbu  waktu  (t)  maka  dikatakan  deret  stasioner dalam mean. Pengujian  stasioner  terhadap mean dapat menggunakan Augmented Dickey-Fuller (ADF) Test.
1 t t Z Z      
Hipotesis : H 0 : 0   (terdapat akar unit, variabel Z tidak stasioner) H 1 : 0   (tidak terdapat akar unit, variabel Z stasioner)
Statistik uji : ˆ ˆ ( ) hitung t SE   
dengan ˆ ( ) SE  merupakan standart  error dari  nilai  taksiran  daerah penolakan : tolak H 0 jika |t hitung | > t tabel.
Apabila  suatu  data  tidak  stasioner  terhadap  varian maka  untuk  mengatasi  ketidakstasioneran  dalam  varians dapat  dilakukan  dengan  transformasi  data.  Secara  umum transformasi Box-Cox untuk beberapa nilai λ disajikan dalam Tabel 1 [2]
Tabel 1. Tabel Transformasi Box-Cox Nilai Estimasi Transformasi -1,0 t Z / 1 -0,5 t Z / 1 0,0 t Z Ln 0,5 t Z 1 t Z

## C. Identifikasi Model VAR
Menentukan  model  VAR  yang  terbentuk  dari  data atau  identifikasi  model multivariate  time  series salah satunya  dapat  dilakukan  dengan  melihat  pola  dari  fungsi Cross-Correlation  Matrix  Function dan Partial  Cross- Correlation Matrix Function.
1. Cross-Correlation Matrix Function
Pada  penerapannya,  bentuk  matriks  dan  grafik  akan semakin  kompleks  karena  meningkatnya  dimensi  vektor, sehingga  untuk  megatasi  hal  tersebut  dapat  digunakan metode koefisien dengan menggunakan tanda (+),(-), dan (.). Tanda  (+)  mengartikan  bahwa  nilai yang lebih besar dari  2  kali  estimasi standard  error yang  menunjukkan hubungan korelasi positif, tanda (-) menotasikan nilai yang  lebih  kecil  -2  kali  estimasi standard  error dan menunjukkan  korelasi  negatif.  Dan  tanda  (.)  menunjukkan yang  terletak  diantara kali  estimasi standard error sehingga tidak adanya hubungan [3].
2. Partial Cross-Correlation Matrix Function
Partial Autocorrelation Matrix Function dapat digunakan untuk  menentukan  orde  model  AR.  Menurut  Tiao  &  Box, partial  autocorrelation  matrix  function pada  lag k yang dapat  dinotasikan  sebagai [3]. Penentuan  orde  AR pada  serangkaian  observasi  akan  disimbolkan  dengan  tanda (+),(-),  dan  (.)  dengan  intepretasi  tanda  yang  sama  seperti yang dijelaskan pada cross correlation matrix function
D. Estimasi Parameter Model
Estimasi  parameter  pada  model  VAR  ( p )  dapat menggunakan  metode Least  Square [4].  Pengujian  pada metode ini menggunakan hipotesis sebagai beriku
H 0 : 0 i   (parameter tidak signifikan)
H 1 : 0 i   (parameter signifikan), di mana i = 1, 2 ,..., p Statistik uji :
ˆ ˆ ( ) i hitung i t SE      i SE  merupakan standart error dari taksiran i  . Tolak  H 0 apabila  |t hitung |  >  t α/2,(n-p) sehingga  diperoleh keputusan parameter signifikan.

## E. Pemeriksaan Model
Apabila telah mendapatkan model yang signifikan maka selanjutnya dilakukan pemeriksaan diagnostik untuk mendapatkan  model  terbaik  dan  beberapa  model  yang parameternya signifikan [4]. hipotesis pengujian: H 0 :   1 2 0 , , , k p R R R  R (residual model VAR ( p ) white noise ) H 1 : 0 k  R (residual model VAR ( p ) tidak white noise )
Statistik uji :
  1 1 0 0 1 k k k k i Q n tr C C C C      
Residual dikatakan white noise apabila nilai
    2 2 k K k p Q   
Pemeriksaan  diagnostik  selanjutnya  adalah  pemeriksaan apakah residual data mengikuti distribusi normal multivariat. Metode  ini  digunakan  untuk  mengetahui  beberapa  variabel yang  memiliki  hubungan  apakah  mengikuti  distribusi normal multivariat atau tidak [5].
H 0 :    residual  karakteristik  kualitas  berdistribusi  normal multivariat
H 1 :  residual  karakteristik  kualitas  tidak  berdistribusi  normal multivariat
Statistik uji :
* 0 1 1 m i W W m    tolak H 0 jika W * < ; , n m c  sehingga dapat dikatakan bahwa residual karakteristik kualitas tidak berdistribusi normal multivariat.

## F. Pemilihan Model Terbaik
Menentukan  model  terbaik  dapat  menggunakan  nilai minimum dari Akaike Informarion Criterion (AIC) [2].
    2 2 ln p pm AIC p n   
di  mana p  merupakan  determinan  dari  matriks  varian kovarian  residual,  m  adalah  jumlah  karakteristik  kualitas, dan n merupakan banyaknya observasi.
G. Dependensi Variabel
Uji dependensi untuk data multivariat dapat menggunakan  uji Bartlett  Sphericity ,  metode  ini  digunakan untuk  mengetahui  apakah  terdapat  korelasi  antara  tiap variabel  dengan  variabel  yang  lain.  Untuk  menguji dependensi  antar  residual  variabel  dapat  menggunakan hipotesis dan pengujian sebagai berikut [6].
Hipotesis : H 0 : R = I (Karakteristik kualitas saling independen)
(2) (3) (4)
(5)
(6)
(7)
H 1 : R  I (Karakteristik kualitas saling dependen)
Daerah kritis: 0 H ditolak jika 2 1 ( ; p(p 1)) 2 u χ     Statistik Uji :   2p 5 u n 1 ln 6       R
Jika gagal tolak H 0 maka  karakteristik kualitas saling independen  sehingga  asumsi  dependensi  antar  variabel  tidak terpenuhi.

## H. Grafik Kendali M
Grafik kedali M digunakan untuk memonitor variabilitas secara  multivariat  dengan  pengamatan  individu  yang didasarkan  pada successuve  difference. Setelah  diperiksa varians dianalisis menggunakan T 2 Hotteling.
I. Grafik Kendali T 2 Hotteling
Grafik  kendali  T 2 Hotteling adalah  Grafik  kendali yang digunakan  untuk  memonitor mean proses  pada  kasus multivariat.  Grafik  kendali  ini  mendeteksi  pergeseran mean proses  dengan  menggunakan  vektor mean sampel  dan matriks kovarians.
J. Grafik Kendali Multivariate Exponentially Weigthted Moving Average (MEWMA)
Grafik  kendali multivariate  exponentially  weighted moving  average (MEWMA)  digunakan  untuk  mendeteksi apakah  ada  perubahan  pada  nilai mean proses.  Salah  satu keunggulan  dari  grafik  kendali  MEWMA  adalah  grafik kendali ini robust terhadap asumsi distribusi normal [7].   1 t t t l       M z M z t merupakan  vektor mean sampel  dan  λ  adalah  pembobot dengan  rentang  nilai  antara 0 1    dan M 0 = 0 . P engamatan dikatakan keluar dari batas kendali atau out of control ketika  nilai 2 T t lebih  besar  dari batas kendali atas. Didapatkan  dari  hasil  simulasi  yang  disesuaikan  dengan besarnya ARL hingga nilainya konvergen [7].

## K. Grafik Kendali Kombinasi MEWMA
Grafik kendali kombinasi MEWMA merupakan pengembangan  dari  grafik  kendali  MEWMA. Grafik kendali ini mengendalikan variabilitas proses dan mean proses dalam satu  grafik  kendali.  Perbedaan  antara  grafik  kendali MEWMA  standart  dan  grafik  kendali  MEWMV  dengan grafik  kendali  kombinasi  MEWMA  dalam  memonitor variabilitas  proses  dan mean proses  adalah  dari  titik pengamatannya.  Titik  mengamatan  pada  grafik  kendali kombinasi MEWMA menggunakan data standardized . Untuk  titik  pengamatan  pada mean proses  fase  I  adalah sebagai berikut:
    1 1 1 2 1 2 , ,..., , ,..., T z z z z z z z k k k kp z k k kp M c E E E E E E     
sedangkan untuk memonitor mean proses pada fase II :     1 1 1 2 0 1 2 , ,..., , ,..., T z z z z z z z k k k kp z k k kp M c E E E E E E   

##  
Untuk  meningkatkan  kemampuan  memonitor  variabilitas proses,  bentuk  kedua  dari  statistika  MEWMA  berdasarkan pada 2 Z ki E disebut  grafik  kendali  M 2 Z 2 untuk  memonitor variabilitas proses pada fase I

##   
       2 2 2 2 2 2 2 2 1 2 1 2 1 1 2 , 2 , ,..., ,..., T z Z Z Z Z Z Z k k k kp z k k kp M n c E E E E E E     
untuk memonitor variabilitas pada fase II
          2 2 2 2 2 2 2 2 1 2 0 1 2 1 1 2
, 2 , ,..., ,..., T z Z Z Z Z Z Z k k k kp z k k kp M n c E E E E E E     
Matriks  varian  kovarian  yang  digunakan  dalam  fase  II merupakan matriks varian kovarian pada fase I yang telah in control [8].

## L. Menentukan Batas Kendali Kombinasi MEWMA
Penentuan  grafik  kendali  kombinasi  MEWMA  dilakukan menggunakan  pendekatan  kriteria Average  Run  Length (ARL).  Pendekatakan  ini  melakukan  simulasi  ARL,  dimana ARL  merupakan  rata-rata  jumlah  titik  yang  harus  diplotkan hingga  mendapatkan  titik out  of control pertama [7]. Untuk mendapatkan  nilai ARL 0 dilakukan  langkah-langkah  sebagai berikut:
1. Membangkitkan  data  sebanyak  10.000  data  dimana parameter yang digunakan menyesuaikan dengan data asli.
2. Mengambil  sampel  random  sebanyak  1000  data  dari 10.000 data bangkitan.
3. Menghitung  nilai  statistik z k M dan  statistik 2 2 Z k M yang dibandingkan dengan batas kendali atas.
4. Menetapkan  batas  kendali  atas  untuk  masing-masing
grafik kendali z k M dan 2 2 Z k M 5. Menentukan Run  Length (RL),  jumlah  titik  pengamatan yang  harus  diplotkan  sampai  titik  pertama  melewati  batas kendali atas.
6. Menghitung ARL 0 yaitu rata-rata RL untuk N iterasi. 7. Simulasi berhenti jika diperoleh nilai ARL 0 = 370

## M. Analisis Kapabilitas Proses
Analisis  kapabilitas  proses  digunakan  untuk  menaksir kemampuan  proses  dalam  memenuhi  spesifikasi  yang  telah ditetapkan.  Suatu  proses  dikatakan  kapabel  apabila  proses dalam keadaan terkendali, berada dalam batas spesifikasi, dan memiliki  nilai  akurasi  dan  presisi  yang  tinggi.  Pada  data multivariat  untuk  memperoleh  nilai  kapabilitas  proses menggunakan persamaan berikut:
1 p MC m i i pi W C   
1 pk MC m i i pki W C   
di mana W i merupakan pembobot dengan 1 1 p i i W    [9].

## III. METODOLOGI PENELITIAN

## A. Sumber Data
Data yang digunakan dalam penelitian ini adalah data sekunder yang diperoleh dari hasil uji pengolahan limbah air PT Japfa  Comfeed  Indonesia  Tbk  Unit  Margomulyo  bagian Health,  Safety,  and  Enviroment  (HSE)  Departmen  HRD  Jl. Margomulyo  No.  36-38  Surabaya.  Data  merupakan  sampel yang  diukur  setiap  hari  Senin  hingga  Sabtu  dengan pengambilan  sampel  yang  dilakukan  setiap  pagi,  data  yang digunakan mulai Januari 2018 – Februari 2019. Sampel yang diukur  adalah  limbah  air  yang  diambil  dari  tank clarifier tahap  setelah  pemberian  larutan  HCl  dan  PAC.  Data  peda penelitian  ini  dibagi  menjadi  dua  fase,  dimana  pada  fase  I merupakan  pengamatan  proses  pengolahan  limbah  pada
(8) (9) (10)
(11)
(12) (24) (25)
bulan  Januari  2018  –  6  Agustus  2018  dengan  total  180 pengamatan  dan  fase  II  merupakan  pengamatan  proses pengolahan  limbah  mulai  7  Agustus  2018  –  28  Februari 2019.

## B. Variabel Penelitian
Variabel yang digunakan dalam penelitian ini dapat dilihat pada  Tabel  3,  dimana  variabel-variabel  tersebut  merupakan karakteristik  kualitas  yang  digunakan  PT  JCI  untuk mengukur limbah air.
Tabel 2. Variabel Penelitian Variabel Keterangan Satuan Spesifikasi Z 1 Total Dissolved Solid (TDS) mg/l Maks 2000 mg/l Z 2 pH - 6-9

## C. Langkah Analisis
Langkah analisis data yang dilakukan dalam penelitian ini adalah sebagai berikut.
1. Cek  autokorelasi  menggunakan  plot  ACF  untuk mengetahui  apakah  terdapat  hubungan  antar  observasi karakteristik kualitas limbah air.
2. Melakukan  pemeriksaan  stasioneritas dalam mean dengan menggunakan ACF dan Augmented Dickey – Fuller Tes, pemeriksaan  stasioner  dalam  varian  menggunakan  jika tidak stasioner dalam varian akan digunakan transformasi Box-Cox pada Tabel 2.1
3. Mengidentifikasi  model  VAR  menggunakan partial autoregression matrices.
4. Estimasi parameter dari model VAR.
5. Pemeriksaan diagnostik residual white noise dan distribusi normal multivariat.
6. Pemilihan  model  terbaik  dengan  memilih  nilai  AIC  yang minimum.
7. Melakukan pengujian dependensi resiudal antar karakteristik kualitas limbah air menggunakan uji Bartlett Sphericity untuk  mengetahui  apakah  terdapat  hubungan antar variabel kualitas limbah.
8. Melakukan  pengendalian  kualitas  menggunakan  metode kombinasi  MEWMA,  memonitor  variabilitas  proses menggunakan  grafik  kendali  M 2 Z 2 dan mean proses menggunakan grafik kendali MZ.
a. Menentukan  batas  kontrol  atas  yang  akan  digunakan pada  grafik  kendali  fase  I  pada masing-masing grafik kendali  M 2 Z 2 dan  MZ  dengan  alogatima  sebagai berikut:
b. Membangkitkan  data  sebanyak  10.000  data  yang memiliki parameter yang sama dengan data asli.
c. Memonitor  variabilitas  proses  pada  fase  I  dengan batas  kontrol  dari  hasil  simulasi  menggunakan  grafik kendali    M 2 Z 2 dan mean proses  menggunakan  grafik kendali MZ.
d. Apabila  ditemukan  pengamatan  yang out  of  control maka  dilakukan  identifikasi  penyebab  pegamatan out of control . Apabila penyebab out of control diketahui maka  pengamatan  yang out  of  control tidak dimasukkan  dalam  perhitungan  lalu  membuat  grafik kendali baru tanpa pegamatan yang out of control .
e. Apabila  pada  fase  I  grafik  kendali  yang  memonitor variabilitas proses dan mean proses sudah berada dalam batas  signifikan  atau  terkontrol  secara  statistik  maka dapat dilanjutkan pada pengendalian kualitas fase II.
f. Parameter in  control pada  fase  I  yaitu  Σ Z0 untuk memonitor  mean  proses  dan 2 Z0  untuk  variabilitas proses  digunakan  pada  fase  II  untuk  megetahui kestabilan proses pada fase II.
9. Melakukan  analisis  kapabilitas  proses  untuk  mengetahui apakah  proses  pengolahan  limbah  air  telah  sesuai dengan spesifikasi yang ditentukan.
10. Menarik kesimpulan dan saran.

## IV. ANALISIS DAN PEMBAHASAN

## A. Karakteristik Data
Karakteristik  data  dari  karakteristik  kualitas  pengolahan limbah air PT JCI adalah sebagai berikut:
Tabel 3. Karakteristik Kualitas Limbah Air PT JCI Karakteristik Kualitas Batas Spesifikasi Rata- rata Varian Min Max TDS Max 2000 1634 73441,6 990 2480 PH 6-9 7,2655 0,1068 6,7 9,2
Berdasarkan  Tabel  3  menunjukkan  bahwa  rata-rata  zat padat terlarut atau TDS dalam limbah air PT JCI Tbk dengan periode  pengolahan  limbah  Januari  2018  –  Februari  2019 berada dalam baku mutu yang ditentukan yaitu dibawah 2000 mg/l.  Penyebab  kandungan  TDS  yang  tinggi  ini  adalah adanya  endapan  yang  berlebih  pada  tangki  sendimentasi sehingga  menyebabkan  zat  terlarut  pada  limbah  meningkat. Karakteristik  kualitas  pH  dapat  dikatakan  bahwa  pH  dari limbah  air  yang  dihasilkan  oleh  PT  JCI  Tbk  netral,  dengan nilai  varian  yang  kecil  yang  artinya  bahwa  data  pH  pada limbah  air  bersifat  homogen.  Penyebab  dari  pH  yang  tinggi adalah karena kadar larutan HCl dan PAC dalam tangki lebih sedikit  dibandingkan  dengan  limbah  yang  masuk  sehingga tidak dapat menetralisir limbah yang masuk.

## B. Pemeriksaan Autokorelasi
Asumsi pertama dalam pengendalian kualitas menggunakan grafik kendali adalah tidak adanya autokorelasi atau  hubungan  antar  observasi,  adanya  kasus  autokorelasi pada  data  dapat  menyebabkan  adanya false  alarm pada grafik kendali. Pengecekan autokorelasi pada data karakteristik  kualitas  limbah  air  dapat  dilihat  dari  plot  ACF masing-masing data karakteristik kualitas.
60 55 50 45 40 35 30 25 20 15 10 5 1 1,0 0,8 0,6 0,4 0,2 0,0 -0,2 -0,4 Lag Au to co rr el at io n Autocorrelation for TDS Gambar 1. Plot ACF TDS
60 55 50 45 40 35 30 25 20 15 10 5 1 1,0 0,8 0,6 0,4 0,2 0,0 -0,2 -0,4 Lag Au to co rr el at io n Autocorrelation for PH

## Gambar 2. Plot ACF pH
Berdasarkan  Gambar  1  dan  Gambar  2  menunjukkan bahwa  beberapa  lag  ACF  TDS  dan  pH  keluar  dari  batas confidance  interval ,  sehingga  dapat  dikatakan  bahwa  pada data  kedua  karakteristik  kualitas  limbah  air  terdapat autokorelasi pada data.

## C. Pemodelan Vector Autoregressive
1. Stasioneritas Varian
Asumsi  dalam  analisis time  series salah  satunya  adalah data  harus  memenuhi  stasioner  terhadap  varian,  untuk mengetahui data telah stasioner terhadap varian dapat dilihat dari  nilai rounded  value pada Box-cox plot,  apabila  data belum stasioner terhadap varian maka dilakukan transformasi sesuai dengan Tabel 1.
Tabel 4.
Stasioneritas Varian
Karakteristik Kualitas Rounded Value Batas Bawah Batas Atas TDS 1 0,57 1,69 pH 1 0,46 1,34
Setelah  dilakukan  transformasi  dapat  dilihat  pada diperoleh  hasil  bahwa  karakteristik  kualitas  pH  telah stasioner  dalam  varian  dengan  nilai rounded  value 1. Sehingga  dapat  dilakukan  analisis  lanjutnya  yaitu  stasioner mean.
2. Stasioneritas Mean Stasioner dalam mean artinya data bergerak disekitar nilai tengah. Salah satu cara untuk mengetahui data telah stasioner terhadap mean atau  tidak  dapat  menggunakan  uji  ADF, apabila  data  tidak  memenuhi  stasioner  dalam mean maka dilakukan differencing pada data. Uji ADF dilakukan setelah data  memenuhi  asumsi  stasioner  terhadap  varian.  Hasil analisis stasioner mean adalah sebagai berikut:
Tabel 5. Stasioneritas Mean Karakteristik Kualitas p-value TDS 0,01069 pH 0,01000 Tabel 6. Nilai AIC Model Time Series AIC Model Time Series AIC VAR[1] 9,0105 VAR[1,3,5] 8,975 VAR[3] 10,698 VAR[1,3,8] 9,021 VAR[5] 11,106 VAR[1,3,10] 8,991 VAR[8] 11,502 VAR[3,5,8] 10,591 VAR[10] 11,659 VAR[3,5,10] 10,607 VAR[20] 11,728 VAR[5,8,10] 10,994 VAR[1,3] 9,001 VAR[1,5,8] 9,003 VAR[1,5] 8,975 VAR[1,5,10] 8,985 VAR[1,8] 9,031 VAR[1,8,10] 9,016 VAR[1,10] 8,996 VAR[3,8,10] 10,617 VAR[3,5] 10,635 VAR[1,3,5,8] 8,991
VAR[3,8] 10,623 VAR[1,3,5,10] 8,977 VAR[3,10] 10,643 VAR[3,5,8,10] 10,592 VAR[5,8] 11,017 VAR[1,5,8,10] 9,001 VAR[5,10] 11,011 VAR[1,3,8,10] 9,008 VAR[8,10] 11,437 VAR[1,3,5,8,10] 8,995
Pada  masing-masing  variabel  diperoleh  nilai  pvalue  lebih kecil dari nilai alpha, sehingga diperoleh kepurusan tolak H 0 . Kedua variabel telah stasioner terhadap mean.
3. Penentuan Orde VAR
Menentukan  orde  VAR  dari  data time  series dapat  dilihat melalui  plot  MPACF  yang  terbentuk.  Data  yang  telah memenuhi asumsi stasioner terhadap varian dan mean. Hasil analisis  menunjukkan  MPACF  memiliki  beberapa  lag terdapat beberapa lag yang signifikan dengan adanya simbol- simbol pada lag ke 1, 3, 5, 8, 10, dan 20.
Penentuan orde VAR dapat dilakukan dengan memilih nilai Akaike  Information  Criterion (AIC)  yang  minimum  pada setiap model yang terbentuk :
Nilai  AIC  pada  masing-masing  kombinasi  model  yang terbentuk dari lag yang signifikan, untuk memperoleh model yang  terbaik  maka  memilih  model  yang  memiliki  nilai  AIC paling rendah dari semua model yang terbentuk. Berdasarkan hasil  analisis  pada  Tabel  4.4  model  yang  memiliki  nilai  AIC paling kecil adalah VAR([1,5]) 4. Estimasi Parameter
Estimasi  parameter  dilakukan  untuk  mengetahui  apakah parameter dalam model yang terbentuk adalah signifikan.
Tabel 7. Estimasi Parameter Model VAR([1,5,8]) Karakteristik Kualitas Parameter Nilai Estimasi P- value Variabel TDS AR1_1_1 0,923 0 TDS(t-1) AR1_1_2 4,765 0,461 PH(t-1) AR5_1_1 0,043 0,264 TDS(t-5) AR5_1_2 3,08 0,64 PH(t-5) AR8_1_1 -0,015 0,651 TDS(t-8) AR8_1_2 7,317 0,263 PH(t-8) pH AR1_2_1 0 0,435 TDS(t-1) AR1_2_2 0,717 0 PH(t-1) AR5_2_1 0 0,16 TDS(t-5) AR5_2_2 0,155 0 PH(t-5) AR8_2_1 0 0,494 TDS(t-8) AR8_2_2 0,019 0,631 PH(t-8)
Tabel  7  menunjukkan  bahwa  terdapat  3  parameter  yang memiliki  nilai  p- value lebih  kecil  dari  alpha  0,05  sehingga menghasilkan  keputusan  tolak  H 0 yang  artinya  terdapat  3 parameter  yang  signifikan  dari  12  parameter  yang  ada. Menurut  Hyndman  dan  Kostenko  pengujian  signifikansi  dan asumsi  dapat  diabaikan,  sebab  yang  diutamanakan  adalah kebaikan model untuk memperoleh ramalan yang tepat [10]. sehingga  model  yang  digunakan  adalah  VAR([1,5,8])  dan residual dari model akan dikendalikan sebagai gambaran dari proses pengolahan limbah yang sebenarnya.
4. Pengujian Asumsi Residual Pemeriksaan autokorelasi dapat dilakukan  dengan pemeriksaan secara visual menggunakan plot ACF. Sehingga plot ACF residual untuk model VAR([1,5,8]) (Gambar 4 dan 5):
Gambar  3  dan  Gambar  4  menunjukkan  masih  terdapat autokorelasi  pada  beberapa  lag,  namun  residual  dari  model VAR([1,5,8])  ini  merupakan  residual  yang paling mendekati
dengan  asumsi white  noise sehingga  residual  dari  model  ini digunakan untuk pengendalian proses pengolahan limbah air.
Pada  analisis  distribusi  normal  multivariat  menunjukkan hasil  dari  analisis  Shapiro-Wilk,  setalah  dianalisis  diperoleh nilai  p- value sebesar  2,2e-16  mendekati  nilai  nol  yang  lebih kecil  dari  nilai  alpha  0,05  sehingga  residual  karakteristik kualitas tidak berdistribusi normal multivariat.
D. Pengendalian Kualitas Limbah Air PT JCI Tbk.
1. Pengendalian Kualitas Limbah Air Menggunakan Grafik
Kendali T 2 Hotteling Fase I 60 55 50 45 40 35 30 25 20 15 10 5 1 0,4 0,2 0,0 -0,2 -0,4 Lag A ut oc or re la ti on Autocorrelation for Resi TDS Gambar 3. Plot ACF Residual TDS. 60 55 50 45 40 35 30 25 20 15 10 5 1 0,4 0,2 0,0 -0,2 -0,4 Lag A ut oc or re la tio n Autocorrelation for Resi pH
Gambar 4. Plot ACF Residual pH
Analisis  ini  digunakan  untuk  memonitor  hasil  dari  proses pengolahan  limbah  air  di  PT  JCI  Tbk  secara  multivariat. Untuk  memonitor  variabiltas  proses  dengan  pengamatan individu menggunakan grafik kendali M dan grafik kendali T 2 Hotteling untuk  memonitor mean proses.  Diperoleh  grafik kendali M sebagai berikut:
Berdasarkan  perhitungan  yang  dilakukan  untuk memperoleh grafik kendali M pada Gambar 5 diperoleh nilai batas kendali atas sebesar 7,37 dan batas kendali bawah 0,05 dapat diketahui bahwa terdapat pengamatan yang keluar dari batas  kendali.  Sebanyak  24  pengamatan  yang  keluar  dari batas kendali atas dan 15 pengamatan yang keluar dari batas kendali  bawah.  Sehingga  pada  fase  I  variabilitas  proses belum stabil.
Gambar 5. Grafik Kendai M Fase I
163 145 127 109 91 73 55 37 19 1 30 25 20 15 10 5 0 Sample T² Median=2,10 UCL=17,18 Grafik Kendali T² Hotteling Gambar 6. Grafik Kendai T 2 Hotteling Fase I.
Grafik kendali T 2 Hotteling yang diperoleh pada Gambar 6 menunjukkan  bahwa  nilai  median  2,10  dengan  batas  kendali atas sebesar 17,18 dan batas kendali bawah nol. Pada grafik kendali  T 2 Hotteling dapat  dilihat  bahwa  terdapat  beberapa pengamatan  yang  keluar  dari  batas  kendali,  sehingga mean proses  pada  proses  pengolahan  limbah  fase  I  belum  stabil. Meskipun  grafik  kendali  ini  lebih  sensitif  dalam  menangkap pengamatan yang keluar dari batas kendali namun banyaknya pengamatan  yang  keluar  dari  batas  kendali  pada  grafik kendali  M  maupun  T 2 Hotteling dapat  disebabkan  karena data melanggar dari asumsi grafik kendali tersebut, yaitu data harus berdistribusi normal multivariat.
2. Penentuan Batas Kendali Grafik Kedali Kombinasi MEWMA
Penentuan batas kendali kombinasi MEWMA menggunakan  pendekatan  kriteria Average  Run  Length (ARL).  ARL  adalah  rata-rata  jumlah  titik  yang  harus diplotkan  hingga  didapatkan  titik out  of  control pertama. Dimana  dilakukan  simulasi  untuk  mendapatkan  nilai  ARL 0 370.  Reynolds and Cho [11] dalam penelitiannya melakukan simulasi Average  Steady  State  Average  Time  To  Signal (Average  SSATS) mengemukakan  bahwa  bobot  optimum yang  dapat  digunakan  adalah  0,11989  dan  0,026. Berdasarkan penelitian Reynolds and Cho (2006) pembobotan 0,11989 memiliki nilai Average Time To Signal lebih  kecil  jika  dibandingkan  dengan  pembobotan  0,026. Sehingga  penelitian  ini  menggunakan  bobot  optimum 0,11989  untuk  memonitor  variabilitas  proses  dan  rata-rata proses.
Alogaritma  simulasi  sesuai  dengan  Langkah  Analisis  poin 9c.  Setelah  melakukan  beberapa  simulasi  dengan  bobot 0,11989  simulasi  dilakukan  hingga  masing-masing  grafik kendali  memperoleh  batas  kendali  dengan  ARL 0 370.  Nilai batas  kendali  atas  untuk  memonitor  variabilitas  proses  atau grafik kendali M 2 Z 2 sebesar 238,523 dan batas kendali untuk memonitor mean proses  atau  grafik  kendali  MZ  adalah 10,35.
3. Pengendalian  Variabilitas  Proses  dan  Mean  Proses Pengolahan Limbah Air PT JCI Tbk Fase I
Analisis ini dilakukan untuk memonitor variabilitas proses dan mean proses pada pengolahan limbah air di PT JCI Tbk, dimana fase I merupakan data pengolahan limbah air mulai 2 Januari  2018  –  6  Agustus  2018.  Hasil  analisis  monitoring variabilitas dan mean proses fase I adalah sebagai berikut:
Gambar 7. Grafik Kendali M 2 Z 2 dan MZ Fase I.
Grafik kendal M 2 Z 2 dan MZ pada Gambar 7 dengan batas kendali bawah untuk masing-masing grafik kendali adalah nol dan  batas  kendali  atas  untuk  masing-masing  grafik  kendali adalah  238,523  dan  10,35.  Pada  grafik  kendali  M 2 Z 2 dapat dilihat  bahwa  terdapat  2  pengamatan  yang  keluar  dari  batas kendali,  pengamatan  yang keluar adalah pengamatan ke-161 dan  166.  Sehingga  variabilitas  proses  pengolahan  limbah  air tanggal 2 Januari – 6 Agustus 2018 belum stabil. Pada grafik kendali  MZ  dapat  dilihat  tidak  terdapat  pengamatan  yang keluar  dari  batas  kendali  sehingga  dapat  dikatakan  bahwa mean proses  telah  stabil.  Penyebab  pengamatan  keluar  dari batas  kendali  adalah karena adanya tumpukan endapan pada penampungan  limbah  dari  boiler  dan  kotornya  filter  karbon aktif,  sehingga  dilakukan  pembersihan  pada  penampung limbah  dan  filter  karbon.  Hasil  pengendalian  limbah  air setelah dilakukan perbaikan sebagai berikut:
Gambar 8. Grafik Kendali M 2 Z 2 dan MZ Fase I Perbaikan I.
Setelah dilakukan perbaikan diperoleh grafik kendali yang tertera pada Gambar 8 dapat dilihat bahwa setelah dilakukan perbaikan masih terdapat pengamatan yang keluar dari batas kendali  yaitu  pengamatan  ke-159  pada  grafik  kendali  M 2 Z 2 sehingga variabilitas proses belum stabil. Pada grafik kendali MZ  tidak  ada  pengamatan  yang  keluar  dari  batas  kendali. Berdasarkan  pencacatan  yang  dilakukan  penyebab  dari pengamatan  keluar  dari  batas  kendali  adalah  karena  filter
pada  pasir  silica  kotor,  sehingga  grafik  kendali  baru  sebagai berikut:
Gambar 9. Grafik Kendali M 2 Z 2 dan MZ Fase I Perbaikan II.
Perbaikan  dilakukan  kembali  dengan  mengeluarkan pengamatan  yang  keluar  dari  batas  kendai.  Sehingga diperoleh  grafik  kendali  baru  pada  Gambar  9  dengan  batas kendali  atas  grafik  kendali  M 2 Z 2 238,523 dan grafik kendali MZ 10,35 tidak terdapat pengamatan yang keluar dari batas kendali. Sehingga dapat dikatakan bahwa proses pengolahan limbah  air  PT  JCI  Tbk  fase  I  telah  stabil  baik  variabilitas proses maupun mean prosesnya.
Sehingga diperoleh parameter matriks varian-kovarian yang in control sebagai berikut
Z0 1 0, 07632
0, 07632 1          
2 Z0 1 0, 005824
0, 005824 1
   
   
Parameter in  control matriks  Σ Z0 digunakan  untuk memonitor mean proses  pada  fase  II.  Matriks 2 Z0  menunjukkan  parameter in  control untuk  variabilitas  proses pada fase II. Setelah diperoleh matriks varian-kovarian yang in  control pada  masing-masing  grafik  kendali  yang digunakan untuk membuat grafik kendali fase II. Pada fase II batas  kendali  yang  digunakan  berasal dari batas kendali fase I.
4. Pengendalian  Variabilitas  Proses  dan  Mean  Proses Pengolahan Limbah Air PT JCI Tbk Fase II
Proses pengolahan limbah air yang digunakan pada fase II mulai  tanggal  7  Agustus  2018  –  28  Februari  2019.  Fase  II dilakukan  untuk  memonitoring  proses  variabilitas  dan mean proses  apakah  terdapat  pergeseran  proses.  Batas  kendali yang  digunakan  pada  fase  II  adalah  batas  kendali  yang dibentuk  dari  data  fase  I.  Parameter  yang  digunakan  pada fase II merupakan parameter in control pada fase I sehingga diperoleh grafik kendali fase II sebagai berikut:
Gambar 10. Grafik Kendali M 2 Z 2 dan MZ Fase II.
Gambar  10  grafik  kendali  M 2 Z 2 pada  fase  II  memiliki batas  kendali  atas  238,523  dan  batas  kendali  bawah  nol menunjukkan  bahwa  tidak  ada  pengamatan  yang  keluar  dari batas  kendali  sehingga  dapat  dikatakan  bahwa  variabilitas proses  telah  stabil.  Grafik  kendali  MZ  diketahui  bahwa terdapat pengamatan yang keluar dari batas kendali, sehingga dapat  dikatakan  bahwa mean proses  pengolahan  limbah  air pada  tanggal  7  Agustus  2018  –  28  Februari  2019  belum stabil.

## E. Analisis Kapabilitas Kualitas Proses Pengolahan Limbah Air PT JCI Tbk.
Analisis ini digunakan untuk menaksir kemampuan proses, apakah  proses  pengolahan  limbah  air  telah  memenuhi spesifikasi  yang  telah  ditentukan  atau  tidak.  Untuk  analisis kapabilitas,  proses  yang  digunakan  merupakan  proses  yang sudah in  control sehingga  menggunakan  pengamatan  pada fase  I  yang  sudah in  control .  Untuk  mengetahui  kapabilitas proses yang digunakan adalah nilai Cp dan Cpk.
Tabel 8. Kapabilitas Proses Univariat Karakteristik Kualitas Cp Cpk TDS 5,62 2,33 pH 3,40 2,73
Tabel  8  menunjukkan nilai Pp dan Ppk untuk kapabilitas proses  secara  univariat,  sehingga  dilakukan  perhitungan kapabilitas proses secara multivariat sebagai berikut:
Tabel 9. Kapabilitas Proses Multivariat Karakteristik Kualitas Bobot MCp MCpk TDS 0,5 4,51 2,53 pH 0,5
Perhitungan  kapabilitas  secara  multivariat  ditunjukkan pada  Tabel  9  dimana  menurut  Raissi  [9]  jumlah  pembobot untuk  keseluruhan  karakteristik  bernilai  satu.  Suatu  proses dikatakan  kapabel  apabila  memiliki  nilai  Cp  ≥  1,33  [7]. Dalam  proses  pengolahan  limbah  di  PT  JCI  Tbk  kedua karakteristik  kualitas  memiliki  keutamaan  yang  sama  dalam proses  pengolahan  limbah  air  sehingga  pembobotan  yang diberikan  pada  dua  karakteristik  kualitas  tersebut  sama. Berdasarkan  perhitungan  yang  telah  dilakukan  secara
multivariat proses pengolahan limbah air memiliki nilai aktual yang tinggi sehingga proses pengolahan limbah telah kapabel.

## V. KESIMPULAN DAN SARAN

## A. Kesimpulan
Berdasarkan  analisis  dan  pembahasan  diatas  dapat diperoleh kesimpulan:
1. Pada karakteristik kualitas TDS terdapat pengamatan yang melebihi  dari  batas  spesifikasi,  penyebab  TDS  tinggi karena  adanya  endapan  pada  tangki  sendimentasi sehingga  menyebabkan  zat  terlarut  pada  limbah meningkat.  Sedangkan  pada  karakateristik  kualitas  pH melebihi  batas  spesifikasi  yang  ditentukan  karena sedikitnya larutan HCl dan PAC dalam tangki pengolahan limbah.
2. Pemeriksaan  autokorelasi  data  menunjukkan  bahwa terdapat  autokorelasi  pada  masing-masing  variabel. Model time series yang terpilih adalah VAR([1,5,8]).
3. Hasil  monitor  pada  fase  I  yaitu  variabilitas  proses  dan mean proses  sudah  stabil  atau  terkendali  secara  statistik setelah  diidentifikasi  penyebab  pengamatan  keluar  dari batas  kontrol.  Pada  fase  II  menggunakan  matriks  varian kovarian in  control diperoleh  hasil  bahwa  variabilitas proses  telah  stabil  namun mean proses  belum  stabil.
Kapabilitas  proses  diperoleh  nilai  MCp  4,51  dan  MCpk 2,53  yang  artinya  proses  pengolahan  limbah  air  telah kapabel.

## B. Saran
Berdasarkan  hasil  analisis  dan  pembahasan  yang  telah dilakukan maka saran yang dapat diberikan kepada PT Japfa Comfeed Indonesia Tbk adalah untuk melakukan pembersihan  yang  terjadwal  pada  masing-masing  elemen pengolahan  limbah.  Perusahaan  diharapkan  melakukan kalibrasi  pada  alat  ukur  TDS  dan  pH  secara  rutin agar hasil yang  diperoleh  dapat  menggambarkan  keadaan  yang sebenernya.

## DAFTAR PUSTAKA
[1]
M.  Reynolds  and  G.  Cho,  “Multivariate  Control  Charts  for Monitoring  the  Mean  Vector  and  Covariance  Matrix,”  J.  Qual. Technol., vol. 38, no. 3, pp. 230–252, 2006.
[2]
Wei,  W.  W.  S.  (2006).  Time  Series  Analysis  Univariate  and Multivariate Methods, Second Edition. New York (United State of America) : Pearson Education. Inc.
[3] Tiao,  Box.  (1981)  Coefficient  Methods.  Dalam  Wei,  W.,  W.,  S. (2006).  Time  Series  Analysis  Univariate  and  Multivariate Methods, Second Edition. New York (United States of America) :
Pearson Education, Inc. [4]
Lutkepohl,  H.  (2005).  New  Introduction  to  Multiple  Time  Series Analysis. Germany : Springer.
[5] Johnson,  R.  A.,  &  Wichern,  D.  W.  (2007)  Applied  Multivariate Statistical Analysis. New Jersey : Prentice Hall.
[6] Morrison, D. F. (1990). Multivariate Statistical Methods (3rd ed). New Jersey, USA: McGraw-Hill, Inc.
[7] Montgomery,  D.  C.  (2013).  Statistical  Quality  Control,  Sixth Edition. United States : John Wiley & Sons, Inc
[8] Reynolds,  M.R  &  Stoumbos,  Z.G.  (2008).  Combinations  of Multivariate  Shewhart  and  MEWMA  Control  Charts  for Monitoring  the  Mean  Vector  and  Covarians  Matrix.  Journal  of Quality Technology, vol.40, pp 381-393
[9] Raissi, S. (2009). Dalam Multivariate Process Capability Indoces on  the  Presence  of  Priority  for  Quality  Characteristics.  Jurnal  of Industrial Engineering International Vol. 5, No. 9.
[10] Hyndman,  R.J.,  &  Kostenko,  A.V.  (2008).  Forecasting  Without Significance Test.
[11] Reynolds,  M.R  &  Cho, G.H. (2006). Multivariate Control Charts for  Monitoring  the  Mean  Vector  and  Covariance  Matrix,  Journal
of Quality Technology, 38(3), 230-252.