Pipeline testing

A lot of developer wonder how to test pipelines - see this Discord thread. The best way we've found is to create fixed, known data that when our transform acts on it, we can make reasonable assertions about what comes out.

Synthetic data

We are operating in the healthcare domain. We have a data set of events at hospitals and we want to turn them into a running total of patients who happen to occupy the hospital on any given day. Our test data should be:

precisely x medical events over D days for y different unique patients, distributed evenly over hospitals {a,b,c} where each patient is admitted on day i mod D and discharged i mod w days later, where i is the unique patient id in range [0,y]. 

Finally, the discharge date is null every z patients because we know we have bad data (urgh).

If we turn this natural language into a (Python) code signature, it looks like this: 

def random_inpatient_data(

    num_rows: int,

    num_patients: int,

    providers: list,

    discharge_null_every=97,

    attendance_date_fn=lambda i: random_timestamp(i, datetime.strptime("2021-11-01", "%Y-%m-%d")),

    start_time_fn=random_time_given,

) -> DataFrame:

    """

    Fake but realistic data for tests.

    :param num_rows: Total number of raws for this data set

    :param num_patients: The number of unique patients in this data set

    :param providers: possible providers from white a synthetic site will be created

    :param discharge_null_every: step size for null discharge dates

    :param attendance_date_fn: A function taking a seed and returning an attendance date

    :param start_time_fn: A function taking a date and a seed and generating a timestamp

    :return: Synthetic Inpatient data.

    """

...

There's some things to note. First, that there are sensible defaults so our tests create the data with only the peculiarities salient to its needs. Secondly, the data is random-ish but the same for any given set of parameters. 

It all sounds a little complicated but the whole thing is some two dozen lines that's called repeatedly across the test suite.

The assertions

There are certain invariants we can assert on. Some are straightforward like the hospital occupancy must be always be equal or greater than zero; or that the number of admissions for a day must always be less or equal to the running total for that day.

Some are a little more involved, for instance the sum of daily occupancy deltas over the time frame is zero (ie, everybody who is admitted is ultimately discharged). Obviously the algorithm must be robust and not count zombie patients who appear to be never discharged - remember that "the discharge date is null every z patients" above?

Another invariant is that we should have a reading for all days in a contiguous block. Remember that the input data is a series of events. If 10 patients are admitted on Monday and 5 are dischared on Friday and noting happens on Tuesday to Thursday, do we still have readings for those dates even though nothing happened? (We should)

Crafting the test data also raised some interesting corner cases that we needed to take back to the business analysts. For example. if a patient is discharged the same day they're admitted, do they show up on that day's occupancy numbers or not? If the discharge date is null what do we do with this patient? 

Conclusion

The use of synthetic data is a powerful tool when building automated regression tests. Tools like Deequ can test data quality but require a curated data sets. This is much harder than it sounds especially when their schemas and semantics change. 

 

Creating code that tests the functionality of your pipeline allows you to refactor your code with confidence. Once you try it, you'll never go back.

My "What data science can learn from software engineering" presentation

Dr Chris Monit and I presented this presentation at the London MLOps meetup last week. TL;DR: maximize your chances of a successful delivery in data science by adopting best practices that the software industry has established.

"Think of MLOps as the process of automating machine learning using DevOps methodologies" - Practical MLOps (O'Reilly)

More ML pipeline debugings

The simpler model is sometimes better. 

We used Meta's Prophet to forecast the flow of patients between care providers. However, we had problems. 

First, the errors were pretty high when back testing. The average error over all care providers was about 20%.

Second, the sum of all care providers in a region was sometimes wildly larger than we'd expect. Curiously, the forecasts for the individual providers in the region were reasonably good.

Thirdly, even if the total numbers of patients flowing between care providers were plausible, it wasn't compatible with the model that forecast the number of patients at the care provider.

A simpler model

We tried a model that used just the historical average number of patients for a care provider. This calculation was over all data on the same day of the week and same month of the year (but ignoring the Covid pandemic when figures were weird).

Example: for any given Monday in January, we simply looked at the average flow between two nodes for all Mondays in every January over all non-Covid years. 

This approach yielded a lower error than Prophet - about 16%. See below for what metric we used to quantify the error. 

Odd figures

When we rendered the Prophet forecasts in the front-end application, they didn't seem too bad most of the time. For instance, we predict that 3027 people would go through the emergency departments in all hospitals in region X tomorrow when today it was 2795. OK, not too shocking.

But if we looked at flows with small numbers, the figures looked crazy. For instance, we predicted the number of people being discharged from all hospitals in a region into mental health units would be 15 tomorrow when this week it had actually averaged about 1.

One of the issues was Prophet itself. Prophet timeseries for small numbers may well have predictions below 0. Obviously, the users did not like us predicting a negative number of patients. But if we simply mapped all negative number to zero, we might get surprises.

Let's take some random numbers:

>>> import numpy as np 
>>> xs = np.random.normal(1, 2, 10)
>>> sum(xs)
9.788255911110234

but:

>>> gt0 = lambda x: 0 if x<0 else x

>>> sum(list(map(gt0, xs)))

12.020263767841016

The forecast for individual care providers was considered good even though Prophet would regularly predict, say, 1 person being sent from an emergency department to mental health facilities when often the figure was actually zero. 

But this 100% error was operationally fine for the clinicians because it was such a small number. It only became an issue when summed over all providers in the region. 

Absent data

We thought using the simpler, historical averages would solve this problem but it didn't. So, we took a closer look at the data. An emergency department would only discharge people into a mental health unit roughly once a week. In the absence of such a discharge, there would be no record. Why should there be? The hospitals only record events that happen not events that don't.

For both Prophet and historical averages, we're not feeding information into our model. But the absence of data is itself information! For instance, the averages of two series are different if you choose to ignore zeros:

>>> np.mean([5,0,0,0,0,3,0,0])

1.0

>>> np.mean([5,3])

4.0

Not a revelation when you think about it but backfilling this information via Spark was a bit gnarly. This StackOverflow answer helped. Once we did this, all figures agreed.

Vertices and Edges

Using historical data, the average flows between nodes necessarily summed to the values at the nodes. This was not true for independent Prophet models.

Imagine we have three care providers a, b and c. 

The number of patients that daily flow through them are A, B and C respectively.

Patients flow from a to c and from b to c.
Let's say the flow from a to c is X and the flow from b to c is Y.

The expected number of patients flowing through c is therefore exactly:

E[C] = E[X+Y] = E[X] + E[Y]  

Compare this to Prophet (or any other model that takes just the historical C to forecast future C). In general:

E[C] != E[A] + E[B]

This is important at the front-end where users are complaining that the predicted flows on their screen don't add up.

Which error metric to use?

I was using root mean squared error to calculate errors whereas a colleague was using  mean absolute error. Which is most appropriate must be decided on a case by case basis. RMSE punishes a data set with a large number of outliers. Whereas, "If being off by 10 is just twice as bad as being off by 5, then MAE is more appropriate". Note that "MAE will never be higher than RMSE because of the way they are calculated" - see this SO answer.

Conclusion

Start with the simplest model you can. It's generally easier to implement, cognitively less demanding and chances are that it will give reasonable results anyway. If the results are not good enough, only then investigate shiny new toys.

MLOps: debugging a pipeline

The domain

Healthcare in England is broken down into about 40 regions. For each region, we want to measure the differences in clinical outcomes conditioned on the ethnic and socioeconomic categories of the patients. To do this, we feed the data for each health region into a Spark GLM.

The problem

Everything was fine with our pipeline for six months before it started to blow up with:

Caused by: org.apache.spark.SparkException: Failed to execute user defined function(GeneralizedLinearRegressionModel$$Lambda$4903/0x0000000101ee9840: (struct<type:tinyint,size:int,indices:array<int>,values:array<double>>, double) => double)

Now, before we do a deep dive, the first thing to note is that we have a robust suite of tests that use synthetic data and they are all passing. 

Secondly, the code that was blowing up was used by five other data sets and they were all working fine in production.

If the code seems OK but one path through the ML pipeline was blowing up in code common to other paths, what does this suggest? Well, if it's not the code, there must be something suspicious about the data, right? The tests use synthetic data so of course they would pass.

The investigation

The first course of action when debugging is to take a good, long stare at the error. This might be obvious but many devs pay insufficient attention to it as it's generally a hundred lines of hard-to-read stack trace. This is like a detective who disregards the crime scene because there's too much evidence to collect. 

Anyway, our murder scene was full of Scala and Python stack traces but if we persevere, we find the line that was triggering the error was a call to Dataframe.collect(). This is generally suspicious but on this occasion, we knew we were dealing with a very small data set so this seemed safe. Indeed there were no OOMEs which is the most common problem with calls to collect(). 

But remember Spark is lazily evaluated. It could be something deeper in the stack that is the root cause. So, navigating a few stack frames previous, we see some one-hot encoding of ethnic groups. Hmm, what can go wrong with one-hot encoding? Well, one potential gotcha is when there is only one category, an exception will be raised.

However, this seemed unlikely. We break down ethnicities into only five groups and there are over a million people in each health region. It would be extraordinarily unlikely if there were a region that only had patients of a single ethnicity. 

Time to look at the data.

Any region with such homogenous patient data probably has very little data to begin with so lets count the number of rows per region. And bingo! there it is: a region called null that has a single (white) patient. This was a recent development in the data being fed into the model which explained why things had worked so well for so long.

The offending row comes from upstream data sets curated by a different department entirely so we're still considering what to do. For now, we could apply a band-aid and filter out any regions called null or better still, any region with fewer than a few thousand patients (as otherwise we're likely to get single cohorts).

One model to rule them?

At the end of the day, the code, the model and the data need to be considered holistically. For instance, which data sets you feed into a model must be evaluated beforehand. 

As an example, we also condition on age bands in this particular GLM model so if we were to feed neonatal or paediatric data into the model it would blow up as all patients would fall into the 0-18 age band. Obvious when you think about it but perhaps surprising if you've inherited somebody else's code.

Tips for effective MLOps

Some miscellaneous tips I've discovered after over a year of being hands-on with a clinical ETL pipeline.

Technology

Set up a local dev environment (Git, Pip, IDE of choice, Python environment etc). Being able to test locally cannot be more important. For instance, I've been using Palantir's Foundry but since it's just a wrapper around Spark, you can have a full, locally run, test suite using PySpark.

If you can't get the data on your laptop (GDPR etc), use simulated data. If your model cannot find a signal in the noise you control, it likely won't when it faces real data. Make these tests automated so it can be part of your CI/CD pipeline.

Keep everything as deterministic as you can. Using commands like SQL's first means you have no guarantees what will come out of the sausage machine. How can your QAs test that?

Prefer to write intermediary data sets to persistent storage and load them again in the next stage over processing everything in memory with a single output. It might make the pipeline take longer to run but it gives you something like an "audit log" of the steps in your model.

Have diagnositic transforms in your pipeline. Typically, they make assertions to ensure self-consistency (eg, all patients' hospital length-of-stay must be non-negative; the aggregated total of overnight patients cannot exceed the total number of beds). These transforms output any violations into data sets that will never be seen by the end user, just the MLOps team. Diagnostic transforms have the advantage over the local test in that they can use real data as they are embedded in the production pipeline.

Process

Testing data is harder than testing code so spare a thought for your QAs. Talk to the them long before you cut a line of code and establish whether they can test what you're building. If not, find out what you need to give them.

Talk to data owners, tell them what you're trying to do and ask them how they'd do it. They often tell you things you'd never know a priori - like which data is, ahem, "less reliable". 

Spend some time (a week?) just playing with the data (EDA) so things become more apparent. It's important to know, for instance, if a column called "region" in A&E data set is the region of the patient or the care provider. Since a patient in A&E may be treated the other side of the country (perhaps they were taken ill while on holiday), the two regions might have nothing to do with one another.

Make sure everybody is singing from the same hymn sheet with a ubiquitous language. Calculating the number of weeks a patient waits couldn't be easier, right? Take the number of days between the start and end date and divide by 7, right? But do we round up or down? And MS SQLServer has its own weird way of doing this calculation.

Be iterative. Don't clean data, build features and then try models in a waterfall manner. 

Data

Don't forget the effect of time on your data (data drift). For instance, care provider networks often get restructured so does hospital trust X in 2020 go by a different name today? Patients move around the country so the same person may be associated with more than one locale. Did more data sources come online over time?

How often is the data updated and how stale is it? These are not the same thing. For instance, the data may be updated daily but its content lags behind by several months.

Were dodgy assumptions made on upstream data sets that you need to be aware of? For instance, removing dead people from the the patients data set might have made sense upstream when dealing with vaccination data but not for historical accident and emergency data.

When unifying data sets, watch out for semantically equivalent values that may be lexically different (eg, if one data set uses an ethnicity of "British white" and another that says "White British", then a union of the two will look like two different ethnicities). This problem becomes harder when one data set uses, say, different age bands. How would you reconcile them?

More MLOps: bug hunting

We have several ML models for which we calculate the odds (the coefficients in a linear regression model) and the risks (the difference in probabilities between the factual and the counterfactual). 

Previously, we saw that for one of our six models, the sign of the odds and risks differed when they should always be the same.

Fortunately, we caught this by putting a transform into our pipeline that takes the output of both the linear regression model and factual/counterfactual probabilities and simply checks them.

Part of our ML pipeline. Some transforms create production data, others diagnostic.

In the graph above, anything with _regression_tests appended is a diagnostic transform while everthing else is production data that will ultimately be seen by the customer.

Recently, a diagnostic test was reporting:

difference in risk (0.14285714285714285) does not agree with the sign of the coefficient (-4182.031458784066) for STP QR1 when p-value = 0.988606800326618

Transform containing pipeline diagnostic data.

Aside: there are other violations in this file but they're not unexpected. For instance, there are 42 health regions in England but the data has some random rubbish from the upstream data capture systems. The last row is due to us expecting at least 10 statistically signifcant outputs from our model when we only found 9. That's fine, we'll turn up the squelch on that one. A red flag would be if the number of insights dropped to zero. This happened when we needed the outcome to be binary but accidentally applied a function of type ℕ => {0,1} twice. Oops.

In a previous post, we found the problem was the "risks" and the "odds" transforms, although using the same input file, were consuming that file at significantly different times. The problem manifested itself when most regions for that one model had different signs. But this single diagnostic we see in the violations above is different. Only one region (QR1) differs. All the others in the data set were fine. 

Looking at the offending region, the first things to note was that we had very little data for it. In fact, of the 42 health regions, it had the least by a large margin.

Number of rows of data per health region for a particular model

So, the next step was to build a Python test harness and using these 127 rows and run the data against production code locally. To run this on my laptop is so much easier as I can, say, create a loop and run the model 20 times to see if there was some non-determinsm happening. Apart from tiny changes in numbers due to numerical instability, I could not recreate the bug.

So, like in the previous post, could my assumption that the data is the same for both risk and odds calculations be flawed? The pipeline was saying that the two descendant data sets were built seconds after the parent.

Both risk and odds data descend from the same parent

Could there be some bizarre AWS eventual consistency issue? This seems not the case as S3 read/write consistency now has strong consistency (see this AWS post from December 2020).

Conclusion

Why just one region in one model has a different sign is still a mystery. Multiple attempts to see the same bug have failed. However, since there is so little data for this region, the p-values are very high (c. 0.99) we disregard any insights anyway. We'll just have to keep an eye out for it happening again.

MLOps debugging: an example

In our ML pipeline, we use generalized linear models to calculate the odds of certain clinical outcomes. We showed this to the client but odds were hard for them to understand. "Can we have probabilities instead?" they asked. 

So, having trained the GLMs, we fit the same data and calculate the average probabilities for each cohort. We then bastardise the data to create the counterfactuals For example, socioeconomic status is of interest so let's make everybody the same level (counterfactual) then make our predictions again and compare with the factual results.

Now, one would imagine that if the corresponding GLM coefficient were positive, then the average probability from the counterfactuals would be less than that of the factual data. This is simple maths, right? You dot product your feature vector with the coefficients and plug it into the sigmoid function. A positive coefficient would lead to a bigger product and so a smaller denominator.

We were mostly seeing this but about 10% of our predicitions violated this rule.

Investigation

A transform compares the two approaches (risks and odds) and with a bit of Pandas:

> import pandas as pd

> df = pd.read_csv("/home/henryp/Documents/CandF/fitted_vs_coeffs.txt", delimiter='\t')

> sub = df[["coefficients", "cohort_risk", "counterfactual_cohort_risk", "p_values"]]

> sub["diff"] = sub["cohort_risk"] - sub["counterfactual_cohort_risk"]

> sub = sub.reindex(sub["diff"].abs().sort_values().index)

> sub[-5:]

    coefficients  cohort_risk  counterfactual_cohort_risk  p_values  feature                       model_id               diff

25     -0.520122     0.239911                    0.184961  0.121023  ethnic_category_one_hot_Other a_and_e_cancer         0.054950

23     -0.561247     0.277648                    0.219477  0.006748  ethnic_category_one_hot_Asian a_and_e_cancer         0.058171

44     -0.292537     0.545455                    0.480899  0.821720  ethnic_category_one_hot_Asian wait_18_ophthalmology  0.064556

50     -0.224494     0.307358                    0.229469  0.480723  ethnic_category_one_hot_Black a_and_e_cancer         0.077889

8      -5.340222     0.636364                    0.551042  0.957624  ethnic_category_one_hot_Asian wait_18_ophthalmology  0.085322

we note that the p-values indicate a lack of statistical significance in the coefficient even if the difference between factual and counterfactual probabilities can be pretty large.

But given the same data should produce the same p-values and coefficients each time (modulo floating point issues), how could this happen?

In Generalized Logistic Regression models, the coefficients don't change much on each run, even when shuffling the data. This is true irrespective of p-values. Let's demonstrate.

Shuffling

The Scala code taken from the official Spark docs for a GLR (here) gives this:

glr.fit(dataset.orderBy(rand())).coefficients.toArray.zip(glr.fit(dataset.orderBy(rand())).coefficients.toArray).map{ case (x, y) => x - y }

res24: Array[Double] = Array(1.231653667943533E-16, -2.220446049250313E-16, 2.220446049250313E-16, 0.0, 1.1102230246251565E-16, -2.220446049250313E-16, -1.6653345369377348E-16, 2.220446049250313E-16, 2.220446049250313E-16, -2.220446049250313E-16)

Running this a few times shows the differences are miniscule even though the p-values are pretty large:

glr.fit(dataset.orderBy(rand())).summary.pValues

res25: Array[Double] = Array(0.989426199114056, 0.32060241580811044, 0.3257943227369877, 0.004575078538306521, 0.5286281628105467, 0.16752945248679119, 0.7118614002322872, 0.5322327097421431, 0.467486325282384, 0.3872259825794293, 0.753249430501097)

Any differences in the coefficients between runs are almost certainly due to floating point arithmetic.

Modification

However, things are very different when we modify the data; even a small perturbation can radically change coefficients with high p-values.

Note that there are only 500 rows in this test data set. But let's fit two models where on average we drop one row at random. Then we see things like this:

def compare2(): Unit = {

  val tol         = 0.002 

  val fitted      = glr.fit(dataset.where(rand() > tol))

  val pValues     = fitted.summary.pValues 

  val stdErrors   = fitted.summary.coefficientStandardErrors

  val coeffs_rand = fitted.coefficients.toArray  

  val other       = glr.fit(dataset.where(rand() > tol))

  val diffs       = other.coefficients.toArray.zip(coeffs_rand).map{ case (x, y) => x - y }

  val diffs_pc    = diffs.zip(coeffs_rand).map{ case (x, y) => (x / y) * 100 } ;  

  val meta        = pValues.zip(diffs).zip(diffs_pc).zip(stdErrors).map{ case(((p, d), pc), s) => (p, d, pc, s) }

  println("%-15s%-15s%-15s%-15s".format("p-value", "difference", "% difference", "std error"))

  meta.sortBy(x => -math.abs(x._1)).foreach { case (p, d, pc, s) => println(s"%-15.3f%-15.3f%-15.3f%-15.3f".format(p, d, pc, s)) }

}

scala> compare2()

p-value        difference     % difference   std error      

0.977          0.051          -224.562       0.793          

0.602          0.101          -24.475        0.792          

0.574          0.043          9.546          0.793          

0.551          -0.035         7.432          0.795          

0.459          0.044          -7.145         0.828          

0.456          0.085          14.590         0.776          

0.283          -0.068         -7.897         0.803          

0.275          0.133          -15.230        0.796          

0.134          -0.067         -5.484         0.811          

0.006          0.119          5.188          0.830        

We need run it only a dozen or so times to see that the coefficients with the largest p-value tend to have the largest percentage discrepancies between the two fits and that the sign can change.

The Culprit

These CLI antics lead me to question the "given the same data" assumption. A closer look at the timestamps of the input data sets revealed they were built a week apart (Palantir's Foundry will by default silently fall back to data sets on a different branch if they have not been built on the current branch. This data may be stale). Well, our data slowly grows over time. So, could these small differences be responsible for big swings in the GLM coefficients? I created in the pipeline a single snapshot of the data from which both the risk and odds values derived and there were zero discrepancies. Sweet.

Model Selection - part 2

We want to predict who will receive a vaccination. Although the outcome is binary (the patient did or did not receive a vaccination), we can build a confusion matrix if for each class, we pretend it's a binary problem. So, if we started with the booked_vaccinated class, our confusion matrix looks like:

true positive:  prediction == booked_vaccinated AND label == booked_vaccinated

false positive: prediction == booked_vaccinated AND label != booked_vaccinated

true negative:  prediction != booked_vaccinated AND label != booked_vaccinated

false negative: prediction != booked_vaccinated AND label == booked_vaccinated

Now we have a confusion matrix, we can look at AUCs to test our model. But, is it the only metric?

AUC or Log-Likelihood

Here are some miscellaneous notes from a data scientist who spent a good chunk of his PhD optimizing models:

1. AUC is good for ascertaining whether the model is well specified and realistic. It can tell you which choice of features are better than others.
2. Log-likelihood tells us how well the model fits the data. If the data is constant and the model structure is the same, the LL can tell us which choice of hyperparameters are better than others.

Note that the hypersurface of log-likelihood for a linear model is always convex. That is, there are no local maxima and they're guaranteed to give a global maxima.

Because log-likelihoods are logs of probabilities, even small relative differences are important. Since delta log likelihood is the ratio between probabilities, differences of a few units can be significant.

Guarding again Regressions

By having a stage in our pipeline that invokes the AUC calculations, we can keep an eye on our models as the code is developed. For instance, I have a Palantir @transform that just prints out a single row like this:

AUC calculations for an ML pipeline

If this changes radically, I know I've done something wrong.

ML Ops

"ML Ops is a set of practices that aims to deploy and maintain machine learning models in production reliably and efficiently" according to Wikipedia. In particular, my job includes "reproducibility and diagnostics of models and predictions" (see the "Goals" section). Despite repeatability being the bedrock of good science, it's surprising how few data scientists appreciate the importance of a pipeline that can reproduce results on demand.

ML Ops is an especially interesting field because, unlike most areas of software and data engineering, the mapping is often not one-to-one. For instance, in a typical banking app, an input (purchase of a financial instrument) leads to just one output (the client now owns a financial instrument) modulo any cross-cutting concerns such as logging. But in ML models, this is not always the case. Including a socioeconomic feature from the data set may change the output from the model in subtle ways. Or it may not at all. In fact, how can you be sure that you really did include it? Are you sure there wasn't a bug?

There seems to be a dearth of good resources for ML Ops, so here are some observations I've made.

Ontologies

Ontologies seem a really nice idea but I've never seen them work. This is partly due to paradoxes that derive from trying to be all things to all men. True life story: the discharge data for some hospitals was 01/01/1900 for about 100 patients. This lead to cumulative length-of-stays to be negative. "OK," says the ontology team. "It's only 100 patients out of 10 million so let's remove them". But the other 99 columns for these 100 patients were fine. So, during reconciliation with external systems, the QA team had a devil of a job trying to find why their numbers did not add up. They had no knowledge of the missing 100 data points that were absolutely fine other than their discharge dates.

Governance

Weekly town halls that encourage questions and maximise the pairs of eyes looking at the data while the data cleaning takes place. The format of this is very much like Prime Minister's Question Time. In addition, a regularly updated blog and threads that all teams can read aid knowledge transfer. What you really don't want is the data engineers working in isolation and not communicating regularly with downstream clients like the data scientists and QAs.

Test data

Fake data can sometimes be a problem. Say, I create some data but want to bias it. I make all people with a socioeconomic score less than X be in one class rather than the other and run a logistic regression model on it. All goes well. Then, one day, it is decided to have the socioeconomic score one-hot encoded to represent poor and everybody else. The threshold for this partition happens to be X. Then my tests start failing with "PerfectSeparationError: Perfect separation detected, results not available" and I don't know immediately why (this SO reply points out that this causes a model to blow up).

Non-determinism

If there are any functions that have non-deterministic results (for example F.first when you aggregate a groupBy in Spark) then you are immediately at odds with your QA team. For example, a health provider has demographic data including socioeconomic status. But this can change over time. They want to reduce patients events to just patients. This means a groupBy on the patient ID. But which socioeconomic indicator do we use in the aggregation? 

Plenty of Data Sets

When pipelining data, it's better to persist more data sets than fewer. This facilitates QA as otherwise it's hard to understand what functions acted on a row. The disadvantage is that a pipeline may take longer to run since we're constantly persisting files to disk. Still, the greater visibility you get is worth it.

There is no tool as far as I know to see the lineage of functions on a single row.

Repeatability

I recommend having non-production transforms that examine the model's output. For instance, I'm currently working on a work flow that looks for insights in health data. I have a transform that counts the number of statistically significant insights that also have a significant effect (how you calculate these thresholds are up to you). The point being, after ever pull request, the pipeline is run and the number of insights are compared to the previous run. If they change significantly, there's a reasonable chance something nasty was introduced to the code base.

Healthy looking (auto-generated) numbers of insights from my models

Without a decent platform, it's hard to have repeatable builds (see here and here for horror stories). Now, I have my issues with Palantir's Foundry but at least it has repeatability built in. Under the covers, it uses plain, old Git for its SCM capabilities. What's even nicer, models can be built on a 'data' branch that corresponds to the Git branch, so you won't overwrite the data created by another branch when you do a full build.