Fastest splitter?

Question

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 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.

Answer 1

PositionIndex[ ] seems to take half of the time of GatherBy[ ]. I use the the test case provided, but evaluate the time with GatherBy[ ], after the comment by MarcoB:

splitGB = GatherBy[tabtest, Last]; // AbsoluteTiming // First
0.0653971

split[list_] := Module[{positions},
  positions = PositionIndex[list[[All, -1]]];
  Table[list[[positions[[i]]]], {i, Length[positions]}]
]

splitN = split[tabtest]; // AbsoluteTiming // First
0.0361735

Ordering is not important:

tab2 = RandomSample[tabtest];
splitM = split[tab2]; // AbsoluteTiming // First
0.0342582

splitByM = GatherBy[tab2, Last]; // AbsoluteTiming // First
0.0616525

splitM == splitByM
True

May be you do not need the partitioned list of values. Using the list of position indexes may be enough, and will avoid doubling the storage required for the data.