Consider some dataset

N1 = 13;
N2 = 130*10^3;
N3 = 2*10^4;
tabtest = 
  Join[Join[RandomReal[{0, 3}, {N1, 4}], Table[{1.}, N1], 2], 
   Join[RandomReal[{0, 3}, {N2, 4}], Table[{2.}, N2], 2], 
   Join[RandomReal[{0, 3}, {N3, 4}], Table[{7.}, N3], 2]];

I need to get a list of lists with the same last element. The obvious one is

SplitBy[tabtest,Last];

However, it is relatively slow, taking 0.15 s on my machine. I made another code:

PhaseSpaceSplitter[phasespace_] := 
 Module[{grouped, lengths, positions, ranges, pso, pdgs},
  pdgs = phasespace[[All, -1]];
  If[Length[Union[pdgs]] != 1,
   (*Grouping consecutive identical elements*)
   grouped = Split[pdgs];
   (*Calculating the lengths of each group*)
   lengths = Length /@ grouped;
   (*Generating the start positions*)
   positions = Accumulate[Prepend[lengths, 0]];
   (*Forming ranges*)
   ranges = 
    Transpose[{Most[positions] + 1, Most[positions] + lengths}];
   pso = Take[phasespace, {#[[1]], #[[2]]}] & /@ ranges
   ,
   pso = {phasespace};
   ];
  pso
  ]

splitted1 = SplitBy[tabtest, Last]; // AbsoluteTiming//First
splitted2 = PhaseSpaceSplitter[tabtest]; // AbsoluteTiming//First
splitted1 == splitted2

0.188216

0.0133451

True

I am okay with how fast the last approach works. However, the speed gain is lost once I apply RandomSample on tabtest:

tabtest=tabtest//RandomSample;
splitted1 = SplitBy[tabtest, Last]; // AbsoluteTiming//First
splitted2 = PhaseSpaceSplitter[tabtest]; // AbsoluteTiming//First
splitted1 == splitted2

0.152886

0.243063

True

How to make an efficient code that makes the list of lists?

P.S. I cannot first make the backward sorting SortBy[tabtest,#[[1]]&] since tabtest is a set of columns of another table, and the ordering is important.