set.seed(231)
data1=gtools::rdirichlet(10,c(1,1,4,4,20,20))
data=t(data1)
k=3
seed=5
index=0

distance<- function(vect1, vect2){
  sqrt(sum((vect1-vect2)^2))
}

centers_function<-function(data, grouping){
  center=matrix(0,length(grouping), dim(data)[2])
  for (i in 1:(length(grouping))){

    if(length(grouping[[i]])==1){
      center[i,]=data[grouping[[i]],]
    }else{
      center[i,]=apply(data[grouping[[i]],],2,mean)
    }
  }
  return(center)
}



if(!is.null(seed)){
  set.seed(seed)
}
centers <- data[sample(1:nrow(data), k), ]


#Step 1.
Dist_e_cent=matrix(0,dim(data)[1],dim(centers)[1])
for (i in 1:(dim(data)[1])){
  for (j in 1:(dim(centers)[1])){
    Dist_e_cent[i,j]=distance(data[i,],centers[j,])
  }
}


Ic12=Dist_IC1_IC2(Dist_e_cent)
Ic12_change=Ic12
Group=Ic12[,1]
grouping<-list()
for(i in 1:(max(Group))){
  grouping[[i]]=which(Group==i)#list with the number of the rows of the data matrix that are in the group [[i]]
}

#Step 2. Update the clusters centers.
centers=centers_function(data, grouping)

#Step 3. Initially, all clusters belong to the live set.
LIVE_SET_original1=c(1:length(grouping))
LIVE_SET_original2=c(1:2)



P1=Step4(data, centers, grouping, LIVE_SET_original1, distance, centers_function, Ic12_change, index)
P2=Step4(data, centers, grouping, LIVE_SET_original2, distance, centers_function, Ic12_change, index)

listfinal=list(c(2,1),c(5,6), c(4,3))



test_that("Step4", {
  expect_equal(dim(P1$centers)[1], k,tolerance=1e-6)
  expect_equal(dim(P1$centers)[2], dim(data)[2],tolerance=1e-6)
  expect_equal(dim(P2$centers)[1], k,tolerance=1e-6)
  expect_equal(dim(P2$centers)[2], dim(data)[2],tolerance=1e-6)

  expect_equal(dim(P1$IC1andIC2)[1], dim(data)[1],tolerance=1e-6)
  expect_equal(dim(P1$IC1andIC2)[2], 2,tolerance=1e-6)
  expect_equal(dim(P2$IC1andIC2)[1], dim(data)[1],tolerance=1e-6)
  expect_equal(dim(P2$IC1andIC2)[2], 2,tolerance=1e-6)

  expect_equal(P1$grouping, listfinal,tolerance=1e-6)
  expect_equal(P2$grouping, listfinal,tolerance=1e-6)

  expect_equal(max(P1$Live_set)<=k, TRUE,tolerance=1e-6)
  expect_equal(max(P2$Live_set)<=k, TRUE,tolerance=1e-6)


  expect_equal(max(P1$no_Change)<=dim(data)[1], TRUE,tolerance=1e-6)
  expect_equal(max(P2$no_Change)<=dim(data)[1], TRUE,tolerance=1e-6)
})