Charcoal, 66 64 bytes

ＷＳ⊞υιυＦＬυＦＬθ«Ｊκι¿⬤⊞ＯＫＶＫＫ¬№_λψ»≔⟦ＬＫＡ⟧ηＦＬυＦＬθ«Ｊκι¤X⊞η⁻ＬＫＡΣη»⎚Ｉ⌈✂η¹

Try it online! Link is to verbose version of code. Takes input as a rectangular list of newline-terminated strings. Explanation:

ＷＳ⊞υιυ

Input the array and write it to the canvas.

ＦＬυＦＬθ«Ｊκι

Loop over every cell.

¿⬤⊞ＯＫＶＫＫ¬№_λψ»

If it is not _ and has no _ neighbour then mark this cell as deleted, so it doesn't show up in a count of cells.

≔⟦ＬＫＡ⟧η

Get the number of remaining cells.

ＦＬυＦＬθ«Ｊκι

Loop over every cell.

¤X⊞η⁻ＬＫＡΣη»

Try to flood fill from that cell, undeleting connected deleted cells, and see how many cells were undeleted.

⎚Ｉ⌈✂η¹

Output the largest orthogonally contiguous area.

Python 3, 145 bytes

I,*g=input(),;i=j=I.find(',')
for _ in I:g+={*sum(([*k]for k in g if{i-j,i}&k),[i+1])*(I[i:i+3]+I[i-j::j+1]>=6*'@')},;i+=1
print(max(map(len,g)))

Try it online!

Takes a comma-separated string as input, with . representing inactive areas and @ representing active areas.

Port of my answer to Max Island Area.

05AB1E, 58 bytes

2Fø€ü3}εε˜ŽqÆSèP]©˜ƶ®gäΔ4F¬ašøí}2Fø€ü3}®*εε˜ŽqÆSèà]˜DÙ0K¢à

Input as a matrix of bits, with X=1 and _=0.

Try it online or verify all test cases.

Explanation:

Step 1: Get a matrix of all surrounded squadrats:

2Fø€ü3}    # Create overlapping 3x3 blocks:
2F    }    #  Loop 2 times:
  ø        #   Zip/transpose; swapping rows/columns
   €       #   Map over each inner list:
    ü3     #    Create overlapping triplets
εε         # Nested map over each 3x3 block:
  ˜        #   Flatten it to a list of 9 values
   ŽqÆ     #   Push compressed integer 13457
      S    #   Convert it to a list of digits: [1,3,4,5,7]
       è   #   Get the values at those indices
        P  #   Take their product to check if all five are 1
           #   (aka an "X" surrounded by four "X")
]          # Close the nested maps

Try just step 1 online for all test cases.

Step 2: Flood-fill that matrix: ©˜ƶ®gäΔ4F¬ašøí}2Fø€ü3}®*εε˜ŽqÆSèà] (see one of the many previous 05AB1E flood-fill answers of mine for an explanation - e.g. this one or this one.

Try just steps 1 and 2 online for all test cases.

Step 3: Get and output the largest island:

˜          # Flatten the matrix to a list
 D         # Duplicate it
  Ù        # Uniquify this copy
   0K      # Remove the (potential) 0
     ¢     # Count for each remaining value how many times it occurs
           # (aka get the size of each individual island)
      à    # Pop and keep the maximum
           # (which is output implicitly as result)

See this 05AB1E tip of mine (section How to compress large integers?) to understand why ŽqÆ is 13457.

PS: I tried un-DRY-ing 2Fø€ü3}εε˜ŽqÆSèP] and 2Fø€ü3}®*εε˜ŽqÆSèà], but unfortunately it's 1 byte longer instead of any shorter due to the minor differences ®* and P/à: try it online.

MATL, 17 bytes

2Y6Z+5=4&1ZIXz&XM

Input is a 2D numerical array with 1/0 for active/inactive.

Try at MATL Online! Or verify all test cases.

How it works

       % Implicit input: 2D numerical array
2Y6    % Push [0 1 0; 1 1 1; 0 1 0]
Z+     % 2D convolution, keeping size
5=     % Equal to 5? Element-wise
4&1ZI  % 'bwlabeln' (label connected components) with 4-neighbourhood.
       % Each connected component, as defined by value 1, is assigned a
       % distinct positive integer
Xz     % 'nonzeros': gives column vector of non-zero elements
&XM    % Second output of 'mode': gives number of occurrences of the mode
       % Implicit display

Uiua 0.17.0-dev.1, 16 bytes ^SBCS

/↥⊜⧻./↧⊂≡⬚0↻A₂¤.

Try on Uiua Pad!

Takes input as a binary array. This code ANDs together the array and the four 2D shifts, then computes the lengths of contiguous regions and takes the largest.